# Review of multiscale geometric decompositions in a remote sensing context

**Mariem Zaouali**,

**Sonia Bouzidi**,

**Ezzeddine Zagrouba**

Université Tunis El Manar, Institut Supérieur de l’Informatique, Research Team SIIVA, Laboratoire LIMTIC, 2 Street Abou Rayhane Bayrouni, Ariana 2080, Tunisia

*J. Electron. Imaging*. 25(6), 061617 (Dec 02, 2016). doi:10.1117/1.JEI.25.6.061617

#### Open Access

**Abstract.**
The quest for optimal representations is considered a challenging goal in the field of image processing. This consists of reducing the processing’s complexity while ensuring an efficient reconstruction. An optimal representation should conserve the properties of the image pertaining to smooth content and contours. The multiscale geometric decompositions (MGD) were designed to reach this finality. They were used in many fields and for different purposes, such as feature extraction, detail enhancing, and change detection. A state-of-art of these decompositions is proposed in this paper. We present their theoretical definitions and how they capture the feature of the objects within an image. An overview table is elaborated where we summarize the methods, the data and the different criteria of assessment used in the studied cases. We are interested, particularly, in the use of MGD in a remote sensing (RS) context. Thus, some examples of their applications on RS images are studied. A discussion is presented based on the analyzed cases.

Since 25 years ago, the definition of the wavelet has led to revolutionary achievements in the development of efficient encoding of piecewise regular signals. Their ability to detect singularities much more efficiently than Fourier transform^{1} has contributed to their remarkable success. The wavelet provides an optimal sparse approximation of a wide set of signals and digitalizes the continuous domain transform through the fast algorithmic implementation.

Nevertheless, this latter fails when dealing with singularities along curves, i.e., $C2$ functions, since reconstructing images with these anisotropic features could produce undesirable artifacts. A wavelet has only diagonal, vertical, and horizontal directions, and this limited directionality sensitivity leads to inaccurate results when dealing with multidimensional data. To overcome this drawback, some forms of directionalities have been introduced, such as the steerable pyramid,^{2} the directional filter banks (DFBs),^{3} and the two-dimensional (2-D) directional wavelet.^{4} A complexifyied version of wavelet was also proposed in Ref. ^{5} to tackle the directionality limitation by adding more precision to $C2$ functions. In fact, complex wavelet offers shift-invariance and yields, at each scale, six directional sub-bands oriented at $+15$ and $\u221215\u2009\u2009deg$, $+45$ and $\u221245\u2009\u2009deg$, $+75$ and $\u221275\u2009\u2009deg$ (indiscriminately by ordinary wavelet).

However, it was proven through several applications that these aforementioned methods cannot represent effectively the image singularities and its anisotropic features. A real transition occurred once Candès and Donoho^{6} introduced the second generation of curvelet (CT) in 2004. Since then, other multiscale geometric decompositions (MGD) were introduced to counteract the CT’s weakness, namely redundancy. In this work, we study some of these MGD and represent the characteristics behind their success. We are also interested in their application in the field of remote sensing (RS).

This field has drawn a great interest with the development of multiple types of sensors, leading to a huge heterogeneous amount of data. To overcome the curse of dimensionality while providing a powerful tool for environment monitoring, planning, and decision-making, several techniques have been employed. Convolutional neural network or the deep learning^{7}^{–}^{9} represents a tendency in this field thanks to its ability to automatically discover relevant contextual features in image categorization problems.^{10} In the same context, sparse representation (SR),^{11}^{–}^{13} total variation,^{14}^{,}^{15} and machine learning^{16}^{,}^{17} techniques, to name a few, are also investigated for different purposes, such as enhancing spatial resolution, generating explanations, and extracting knowledge from the images. MGD, the scope of our paper, have been extensively used in several domains, namely pattern recognition^{18}^{,}^{19} and computer vision.^{20}^{,}^{21} Their use was also extended to the RS field, where we need to locate edges of roads, building, rivers, and forest to detect a potential change. Adding to that, MGD provide an analytical treatment of a scene by decomposing it in high frequencies and low frequencies. This essentially helps in fusing information from different sensors and in revealing hidden characteristics indiscriminately using only one sensor image. Therefore, we aim, in this work, at drawing attention to the MGD importance and contributions in the RS field. To the authors’ knowledge, several reviews were elaborated to describe MGD, but few of them focused on their particular use in RS.

This paper is organized as follows: we present, in the first section, an overview of MGD. Then we describe their characteristics and how they have been used in the field of RS. We conclude this paper with a discussion and a conclusion.

The notion of scales was introduced before wavelet emergence to better localize objects in the image observation. Methods like windowing or scale introduction in Fourier transform,^{22}^{,}^{23} Laplacian pyramid (LP), and derivatives of Gaussian (DoG) were the cornerstones of ideas a few years later. The definition of the continuous wavelet was a kind of generalization of previous works in scale incorporation since the newborn transform has offered a localization in time and frequency.

For numerical computation needs, the wavelet was discretized using multiresolution analysis (MRA). The MRA consists of projecting a signal successively onto subspaces $Vj$ (as $j\u2208Z$ increases, the approximation becomes coarser) and yielding both approximation and detail information (detail information is calculated from two successive approximations). Otherwise, MRA stipulates that any signal, in our case an image, could be constructed iteratively by exhibiting different characteristics in every scale. Figure 1 shows how smooth regions are emphasized in finest scale while contours are more likely to be salient in the coarsest ones. Taking advantage of the MRA, wavelet transform can address point-like singularities and offer less redundant information in scales, compared to LP and DoG. Besides, wavelet is characterized by being separable, which means that its 2-D function atoms could be written as the product of two other 1-D functions. Precisely, a wavelet’s atom $\psi (x)=2\u2212j\psi k(2\u2212jx\u2212n)$ is obtained by three tensor products 2-D wavelets as follows:

where $x=(x1,x2)\u2208R2$ refers to spatial coordinates of a given point $x$, $2\u2212j$ refers to dyadic scalings, $j\u2208Z$ refers to scale, $n\u2208N$ is translation factor, $\phi $ is 1-D orthogonal scaling (acts like low-pass filter), $\psi $ is wavelet function (acts like band-pass filter), and $k$ stands for the only three directions supported by wavelets (horizontal, vertical, and diagonal).

The MGD were introduced to further improve the aspect of directionality limitation, whether by using the wavelet itself combined with other treatments or by defining rules for wavelet functions. These methods can be divided into two different families:

- adaptive family, which follows a Lagrangian approach, such as bandlet and
- nonadaptive family, which follows the Eurelian approach, such as CT.

The first family proposes to build an adapted structure of the data. For example, in order to locate singularities, the bandlet decomposes the image using wavelets combined with dyadic decompositions. Each square of the resulted segmentation contains one direction, i.e., a unique piecewise singularity. However, the second family tries to ameliorate wavelet by defining additive mathematical rules on the wavelet functions. This brings more accuracy to detect image content. A debate was settled on the strength of those families and their abilities in image understanding: smooth regions, contours, and texture. According to Ref. ^{24}, the best representation is left an open question since it is data dependent.

In this section, we represent some of the MGD and describe their main characteristics and differences. Particularly, in this paper, we are interested in some MGD that are applied in the RS field. Figure 2 illustrates a map of MGD elaborated based on Ref. ^{24}. In fact, we have written the names of methods that will be discussed in bold, so that the reader can have a clear vision of what is treated in the coming sub/sections.

The first generation of CT was proposed by Donoho and Duncan.^{25} By that definition, this transform extended the ridgelets (built on the ground of 1-D wavelets to detect segment singularity)^{26} to switch from a global scale analysis to fine scale analysis. The idea was to decompose the image into a set of subimages and perform ridgelet transform on each one of them.

This CT transform has not shown great results because of the failure of ridgelets in diagnosing the $C2$ edges. The second generation of CT is then proposed and was less redundant and faster compared to its first generation.^{27} Given a function $f\u2208L2(R2)$, the coefficients $C(j,l,k)$ of the continuous CT transform are defined by Eq. (1):

CT’s atoms are defined by a combination of a translation $xk(j,l)$ and a rotation $R\theta l$ of an angle $\theta l$ of an atom $\varphi j,l,k(x)$:

The atom of CT $\varphi j$ is defined in the frequency domain by the mean of its Fourier transform $\varphi j(\omega )^=Uj(r,\theta )$, which is written in polar coordinates as follows:

To discretize the CT, the authors of this transform propose to switch from rotated tiling to angular tiling and from concentric circles to concentric squares, as illustrated in Fig. 3.

Contourlet is considered as the extension of Candes and Donoho’s work to improve the isotropic criterion of wavelet representation.^{28} This transform aims at obtaining the same frequency space tiling as the CT, without the need to move from continuous to discrete domain. To do so, the authors propose a nonseparable decomposition scheme by applying a DFB, combined with LP [as shown in Fig. 4(a)] and they also propose to decrease the high redundancy information in its sub-bands.

The DFB is designed to capture anisotropic features in the high frequencies of the image by allowing different numbers of direction at each scale while the low frequencies are processed by LP. In fact, as presented in Ref. ^{28}, LP is defined based on orthogonal filters and downsampling by 2 in each iteration, similar to wavelet. The low-pass synthesis filter G, used for LP processing, defines a unique scaling function $\phi (t)\u2208L2(R2)$, $t\u2208Z$ that satisfies Eq. (5):

Otherwise, the DFB offers localization and directionality properties due to the family, defined in Eq. (6):

^{29}discards the sampling step and is characterized by being shift-invariant and by preserving multiresolution criteria.

Shearlets^{30} provide a rich mathematical structure and are optimally sparse in representing the edges within an image thanks to the fact that they form a tight Parseval frame at various scales and directions. They are distinguished from CTs by being directly constructed in discrete domain, which gives them the ability to provide a more efficient multiresolution representation of the geometry.^{31} Besides, shearlet outperforms contourlet by substituting the directional filter with shear filter, which helps in breaking the limitation of directionalities. In the continuous domain, in dimension 2, a shearlet represents an affine system with composite dilations in Eq. (7):

^{32}, the discrete collection of shearlet transform could be defined by two window functions localized on a pair of trapezoids, as illustrated in Fig. 5(a).

In fact, for a given point in the frequency domain $\xi =(\xi 1,\xi 2)\u2208R2^$, $j\u22650$, and $l=\u22122j,\u2026,2j\u22121$, let

^{32}, and $\chi D0$ and $\chi D1$ are, respectively, truncated vertical and horizontal cones. The authors of this transform have also proposed a fast algorithmic implementation and a filter bank architecture which is similar to the one of contourlet transform as shown in Fig. 5(c)].

Wedgelet is proposed by Donoho et al.^{34} This transform decomposes the image iteratively into piecewise constant functions (class of horizon functions). The first step consists of a dyadic recursive decomposition. In each quadtree leaf, wedgelet will search for an “edgel” in order to forward the leaf dyadic decomposition. Figure 6(a) represents a wedgelet decomposition of “cameraman” image and its dyadic decomposition.

The intuitive evolution of this transform is to switch the constant horizon functions by polynomial functions in order to catch more irregular singularities. In fact, this was proposed with the surflet.^{36}Figure 6(b) illustrates how the dyadic decomposition is influenced when polynomial functions are used instead of piecewise constant functions.

The first generation of the Bandlet transform was proposed by Le Pennec and Mallat.^{37} The main idea was to explore the wavelet transform and overcome its failure in detecting anisotropic regularities, such as $C2$ curvatures. According to Ref. ^{37}, Bandlet first performs the grassfire algorithm^{38} to detect edge regions and then they perform a deformation to let it fit wavelet orientations (horizontal, vertical, or diagonal). This version was redundant and failed to ensure good results, mainly in compression applications.^{35}

Thus, the second generation of bandlet is proposed. Respecting the same adaptive approach, this version is considered as an anisotropic wavelet wrapped along the geometry flows.^{39} This version differs from its previous version by being based on quadtree segmentation algorithm. The idea is first to estimate irregularities in the image as a vector field, generated using the edge polynomial function estimation. A wavelet transform is then applied, after that a dyadic segmentation is performed in wavelet domain. This way, wavelet coefficients are calculated along the optimally found direction in each square. The bandlet could also be seen as a polynomial approximation in an orthogonal Alpert basis. Alpert transform is simply a wrapped wavelet transform adapted to a given irregular form, as illustrated in Fig. 7. The coefficients of the bandlet are calculated by following Eq. (8):

There is no doubt that wavelet was a genuine processing tool that has led to major advances in natural image representation and understanding. Since then, many transformations have been defined to overcome the shortage of wavelet capacities by adding directions, shift invariance criteria, and so on. In Table 1, we summarized the characteristics of different transforms mentioned in this review. We have separated the MGD into two families. Their differences reside on the fact that the adaptive family constructs a special data decomposition that can fit each dataset rather than using a predefined system like the nonadaptive families. However, this adaptive form is very intensive in terms of numerical computation.^{1}

Nonadaptive family | |||
---|---|---|---|

Name | Authors | Description | Conception |

CT^{6}^{,}^{25} | Emmanuel Candes and David Donoho | Can characterize $C2$ curvature | Evolution of ridgelet transform |

Cannot characterize oscillating textures directionality sensitive | Redundant | ||

Contourlet^{28} | Minh N. Do and Martin Vetterli | Evolution of isotropic wavelet transform | Use of nonseparable filters |

Can characterize $C2$ curvature | |||

Cannot characterize oscillating textures directionality sensitive sampling discards low frequencies information | Filter banks-based transform | ||

Shearlet^{30} | Demetrio Labate, Wang-Q Lim, and Glenn Easley | Evolution of CT more directionalities than the CT cannot characterize oscillating textures | Substituting the rotation and anisotropic stretch with anisotropic shears |

Adaptive family | |||

Wedgelet^{34} | David Donoho | Use a wavelet quadtree segmentation and grassfire algorithm to detect geometric singularities | A leaf of quadtree supports only horizon function |

Bandlet^{37}^{,}^{41} | Peyre, Le Pennec and Mallat | Extends the isotropic wavelets quadtree segmentation | Wrapping wavelet coefficients along geometric flows |

The second generation is less redundant than the first generation |

But, in general, all these transforms exhibit interesting parameters, which helps to bring attention to different details within an image.

**Scale**: MGD locates, in each scale, a specific feature (edge in the coarsest scale, smooth content in the finest). Keeping some scales while withdrawing others could help remove undesirable information. Indeed, discarding disturbance that could be part of the finest scales helps in providing a better representation of the image**Direction**: Likewise, directions reveal features in the image. We can detect, thanks to this parameter, the main axes by locating the accumulation of strong coefficients magnitude. Consequently, we can enhance some image structures. Figure 8 illustrates a rendering of the content of scales and supported directions of some MGD.**Magnitude**: In transformation domain, the greatest magnitude values could be interpreted as an influent part in the processed image. This is why the magnitude is useful for feature extraction applications. In fact, since the MGD allow the zooming in and out, the content of the image, especially texture, will not be affected by the size of the neighborhood. Several first-order and second-order statistics could be calculated: mean, standard deviation (SD), energy, entropy, contrast, sum of mean, variance, and so on. Consequently, several descriptors could be extracted and could be dimensionally reduced without affecting the discriminative power.^{42}Another interesting use of magnitude is detecting spatiotemporal directions as in Ref.^{43}.

**Frames**: CT, shearlet, and contourlet are a set of functions forming a tight frame. This means that for a given sequence $(\varphi i)i\u2208I$ in Hilbert space $H$, where $I$ is a countable indexing set, $i$ is an element of $I$, and for all $x$ in $H$, there are two constants $0<A\u2264B$:Display Formula $A\Vert x\Vert 2\u2264\u2211i\u2208I|\u27e8x,\varphi i\u27e9|2\u2264B\Vert x\Vert 2.$(9)When $A=B$, this is called a tight-frame, and if $A=B=1$ this is called a Parseval frame. The advantage of having such a function system is to be able to represent a signal as a linear combination of the vectors within the frame in several ways and to reconstruct it by the use of inner products [Eq. (10)]. Let $S=T*T$, where $T$ is the analysis operator and $T*$ is the synthesis operator of a given frame. $T$ and $T*$ will be explained in “sparsity” paragraph:Display Formula $S:H\u2192Hx\u21a6\u2211i\u2208I\u27e8x,\varphi i\u27e9\varphi i.$(10)

The authors of bandlet have also proposed its tight frame version, which is called grouplet.^{44} For further details on this frame version, we refer the interested reader to Ref. ^{45}.

**Redundancy/overcompleteness**: Although the orthonormal/orthogonal bases has interesting properties, they have several weaknesses. They are translation and rotation sensitive, especially when it comes to dealing with multidimensional data. That is why the MGD, such as curvelets, shearlets, and contourlets, are redundant. This means that the vectors of their basis are linearly dependent and exceeded the number of their space dimension. Thus, a given vector in this kind of space would have an infinite set of representations. In fact, the authors of these transforms want to move “cautiously” to overcompleteness without losing the interesting features of the orthonormal/orthogonal basis.**Sparsity**: The aforementioned MGD have two main operators: the analysis $T$ and the synthesis $T*$. Given a signal $x$, MGD analyze $x$ by decomposing it in a set of coefficients:Display Formula $T:H\u2192\u21132(I)x\u21a6(\u27e8x,\varphi i\u27e9)i\u2208I.$(11)The dual operator, the synthesis, helps in recovering a given $x$:Display Formula $T*:\u21132(I)\u2192H[(ci)i\u2208I]\u21a6\u2211i\u2208Ici\varphi i.$(12)Curvelets, contourlets and shearlets use fast algorithms to resolve Eq. (11). A suitable choice of a representation system enables a SR which means that a given signal $x$ could be expressed in a low-dimensional space. This representation is widely used in several approaches such as compressed sensing.^{46}Figure 9 shows the huge difference between representing an image using histograms of a gray level and representing it sparsely using curvelets (as an example of the effect of the use of MGD).

RS images represent an efficient tool to analyze climate change and dynamics of land cover. In the last decades, a great number of new satellites has been launched, allowing a tremendous data availability with improved spatial and spectral resolutions. This has helped in enhancing our understanding and control of our surroundings. Throughout the various RS applications (image fusion, enhancement, super resolution, and so on), we are capable of monitoring the earth’s surface, predicting changes, and preventing disasters.

Whereas, the management of RS images represents a challenging task. As the amount of data is incredibly growing, it is getting more complex to extract knowledge from this type of images. RS data are not only different but also have a rapid velocity and generally need several complicated corrections. That is why researchers seek new approaches to replace the traditional data processing algorithms.

MGD were used in an attempt to solve some of the aforementioned problems. In fact, being able to represent the datasets in a sparse way while describing accurately, the objects in the image make the MGD an inspirational tool to be exploited in RS field. In this section, we are interested especially in the use of wavelet, CT, contourlet, shearlet, and bandlet.

The classification aims at recognizing the class label of a given study area with the aid of ground truth data. Conventional classification algorithms exploit essentially the spectral information of the images.^{48} Nevertheless, it has been proven that incorporating spatial information in the classification process, i.e., taking the contextual information into account, helps significantly in boosting the obtained accuracies. The spectral–spatial combination is possible using MGD, which not only provide a frequency representation but also help in establishing correlations between neighbors in the spatial domain.

- Decomposing the whole image and using its frequency details for classification
- Combining the MGD frequency details with other features obtained using other methods.

We have studied the Ref. ^{49} where wavelet has been used for a feature extraction purpose. Combined with fuzzy hybridization, the wavelet coefficients are exploited for classification. The three multispectral (MS) RS images, two from Indian remote sensing and one from SPOT, are decomposed band-by-band. This technique has essentially helped in accounting for contextual information. The overall classification accuracy ($Mean=79.15%$) using ground truth of the three images has shown that the biorthogonal wavelet exceeded the other wavelet functions.

In Ref. ^{50}, the authors proposed a CT-based approach to extract features from SAR images that can effectively identify the dynamic ice from the consolidated one. The representation of the dynamic ice contains curves in different locations with different widths and lengths and is considered larger than those of consolidated ice representation. The proposed approach consists of extracting patches from the transform domain sub-bands according to a specific size. This later increases according to the CT’s scale in order to capture more information about the curves. Once the size of patches is computed, a feature vector is calculated by Eq. (13):

Using SVM along with different features, the obtained experimentation results proved that the CT-based feature extraction is effective for classification since the dynamic ice area is more accurately classified. In fact, thanks to the probability density function (PDF), we can see clearly that the component representing the different types of ice could be distinguished and thus separated (Fig. 10).

**F10 :**

Estimated PDF for different features: (a) the PDF of SD, (b) gray-level co-occurrence matrix, and (c) CT-based method. The extracted features overlap a lot, which means they will be less effective in classification except for the PDF of CT-based method, which shows that the different types of ices could be distinguished in classification.^{50}

The contourlet is used in Ref. ^{48}, where a performance comparison between wavelet and contourlet is discussed. Based on the fact that wavelet suffers from its shortage of directionality and the fact that contourlet provides directions only in high frequency coefficients, a wavelet-based contourlet transform (WBCT) is proposed and is applied on linear imaging self-scanner (LISS) II, III, and IV. Figure 11 illustrates the proposed transform.

After decomposing the image using the WBCT, wavelet, and contourlet, PCA is applied to the obtained features in order to reduce their dimensionality while removing redundancy and preserving the most discriminant ones. Then a mean vector and a covariance matrix are calculated. Finally, the Gaussian kernel fuzzy C-means classifier is applied, and the obtained overall classification accuracy proves that the proposed decomposition (overall accuracy of $LISS\u2009IV=89.57%$) is better than the wavelet-based (overall accuracy of $LISS\u2009IV=87.62%$) and contourlet-based feature extraction methods (overall accuracy of $LISS\u2009IV=88.38%$).

As we mentioned earlier, the second category combines MGD sub-bands as features describing edges or textures, with features obtained using other methods. For example, the authors of Ref. ^{51} suggest to use CT jointly with morphological component analysis (MCA) and to improve RS classification using hyperspectral AVIRIS and AirSAR images. The idea is about separating a given image into two components (14):

^{52}

This approach was extended in Ref. ^{53}, where the authors combined several methods such as Gabor wavelet and horizontal filters to construct an MCA kernel for feature extraction. The use of such a composite kernel has given better results (overall accuracy of AVIRIS $Indian pines=93.54%$) in characterizing image content compared to minimum noise fraction (MNF) components and helped in enhancing the characterization of the image’s content. This is explained by the fact that the proposed approach combines several methods, such as wavelet and CT.

In a similar fashion, the authors of Ref. ^{54} proposed to combine multiple features pertaining to spectral, texture, and shape, and proposed a multiple feature combining (MFC) framework, as shown in Fig. 12. The spectral feature of a given pixel is elaborated by arranging its digital number in all of the $l$ bands. The texture feature is obtained by applying 2-D Gabor wavelet filter and the shape feature is constructed due by the pixel shape index method.^{55} To calculate the feature vectors using the MFC framework, a single feature-based dimensionality reduction technique is conducted in order to generate the alignment matrix. Then a Lagrangian function is calculated in order to determine the low-dimensional feature space $Y$ (16):

**F12 :**

Multiple features of the airborne data over the mall in Washington DC dataset. (a–c) Spectral feature images in band 36, 52, and 65. (d–f) Gabor texture feature images, with $direction=1$ and $scale=1$, 3, and 5, respectively. (g–i) Shape feature images in direction 1, direction 8, and direction 16, respectively.^{54}

The detection of a change in land cover or in the Earth’s surface is considered as one of the most important applications in RS. In fact, it helps in disaster management, vegetation development, deforestation detection, and urban growth tracking. This application needs a set of multitemporal satellite images. Throughout the studied cases in this section, MGD were extensively used to ensure an accurate detection of affected regions in multitemporal images by enhancing the image’s details, especially the difference images (DIs), where noise information can easily be interfered.

In Ref. ^{56}, the contourlet is used to denoise the images of interest in order to enhance the change detection process. First, it is applied to each single temporal SAR image to preserve its features and edges. The authors also proposed to reduce speckle noise by performing hard thresholding on its high frequency sub-bands using Eq. (17):

After that, markers are extracted from the contourlet domain in order to find the potential course rivers and eliminate the false alarms (content detected as river courses while it is not). After being processed, the SAR images become smoother than the ones before applying the noise and the speckle reduction. In contrast, the edges are preserved and become more visible. That is why the authors affirmed that this contourlet-based method can achieve higher accuracy of river courses detection.

In Ref. ^{57}, the CT is used jointly with MCA and is applied to a high-resolution airbone SAR image to detect changes in urban areas acquisition. The MCA decomposes the image into $K$ different components, as specified in Eq. (18):

To suppress the undesirable information, the authors applied a soft-threshold on the CT and wavelet coefficients. Then they calculated the change map, as illustrated by Eq. (20):

In Ref. ^{58}, the authors detected anomaly in hyperspectral images by, first, using the shearlet to decompose the images into several directional sub-bands at multiple scales. Then, in each sub-band, the background signal is reduced and the fourth-order partial differential equation is applied to brighten up the anomaly. Experimental results with HYDICE HI data show that the presented algorithm can suppress the background, detect the anomaly signal effectively, and outperform the original RX algorithm.^{59}

Another approach is proposed using a detail-enhancing approach using NSCT Ref. ^{60}, where the authors suggested to detect change from multitemporal optical images using a detail-enhancing approach. To do so, they studied and analyzed two change detection methods. The first^{61} is based on combining PCA and K-means, which was efficient in terms of computational time, but since the PCA does not consider multiscale processing, when fusing data, this has yielded false detections. The second method was proposed by Ref. ^{62} and introduces a multiscale framework for change detection in multitemporal SAR images. This method decomposes a DI into $S$ scales by the undecimated discrete wavelet transform (DWT). Based on this decomposition, the authors extracted a feature vector by first computing the intrascale features (sampling local neighborhood using two methods) in each sub-band in a given scale, and then by computing the interscale features, which regroup all the intrascale vectors together.

The same authors improved this work in Ref. ^{63} by proposing the use of CT-DWT (complex wavelet) to capture further directions. The authors also used the Bayesian inference to calculate a threshold to decide whether the details to fuse correspond to changed or unchanged pixels.

In Ref. ^{60}, both works Refs. ^{61} and ^{62}, were discussed. The authors estimated that the use of K-means algorithm, in Ref. ^{61}, can stick to local optima, which can produce false detections, especially when choosing an inappropriate initial centroid. Moreover, they also affirmed that the use of wavelet, in Ref. ^{62}, failed to characterize the geometric details of the DI efficiently. To overcome this weakness, instead of wavelet, the authors used the NSCT to decompose the DI into $S$ scales. Each scale yielded one low pass approximation sub-band $DSA$ and $L$ high-pass directional sub-bands $D1,1H\u2026DS,LH$. The details are extracted from directional sub-bands. The high-pass directional sub-bands serve to extract the details and are thresholded to decrease the amount of noise. The intrascale fusion is performed using max rule, as shown in Eq. (23):

Then an enhanced DI is calculated [Eq. (25)]:

^{64}is performed to calculate the change map. Compared to approaches based on EM, on Bayesian, on PCA, or on multiscale, the proposed approach gives better results since it conserves the geometric details due to the NSCT. The results of experimentations are illustrated in Fig. 14.

**F14 :**

Qualitative change detection results by using different change detection methods. (a) Ground truth change detection map created by manually analyzing the corresponding input images, (b) EM-based approach, (c) Bayesian-based approach, (d) PCA-based approach, (e) multiscale-based approach, and (f) proposed method.^{60}

In Ref. ^{65}, the proposed method applies wavelet on MS imagery in an anisotropic diffusion aggregation. The proposed approach is composed of three steps. We only cite, in depth, the steps where wavelet is used:

- A band selection is first conducted to choose where the texture is best captured. Then, 2-D DWT for multiscale-multidirectional texture extraction is applied.
- Textural and spectral segmentation are performed by anisotropic diffusion in order to reduce noise without blurring inter-region edges as well as creating the desired multiscale low-level primitives.
- Change detection is then performed in two steps:
- classification of the image into forest/nonforest thematic maps through thresholding
- comparison of the classification map, obtained from the previous step, using logical modeling in order to identify changes

The experiments are conducted on Landsat TM and Landsat ETM+ datasets dated 1986 and 2001, respectively. The criteria of assessment used in this work are the transformed divergence measure which statistically determines the adequate wavelet levels to keep, overall accuracy and kappa.

In Ref. ^{66}, a CT-based change detection algorithm is proposed between two co-registered SAR images for natural disaster mapping. After applying the CT on the two images, the coefficients are weighted to suppress noise-like structures using the mathematical relation [Eq. (26)]:

Image fusion represents one of the most important RS applications. It aims at combining two or more images to create another one with enhanced features. This increases the possibility of taking advantage of multisensors images and opens the door for several uses. Generally speaking, fusion can be categorized into four families based on when the fusion rule is applied: fusion at signal level, at pixel level, feature level, and decision level. The MGD belong to the feature level family. In fact, in order to conduct a fusion process with these latter decompositions, we should first convert the intrinsic image properties and get the frequency coefficients. Fusion techniques, based on MGD, suggest several injection models but tend basically to extract the spatial detail information from high spatial resolution images to inject them in the low spatial resolution ones. Compared to the largely used methods like intensity hue saturation (IHS) and PCA, MGD enhance the spatial characteristics of fused images. Thus, they look very clear, have sharp edges, and are free of spectral distortion.^{67} MGD-based fusion are applied abundantly in pansharpening techniques, multisensor fusion, and spatial enhancement.

In Ref. ^{68}, the authors propose a technique to fuse MS satellite images (high spectral and low spatial resolution) with a panchromatic (PAN) satellite image with low spectral and high spatial resolutions. They introduced an improved method of image fusion based on the ARSIS concept using the CT transform. The ARSIS is a French abbreviation for enhancing spatial resolution by injecting structures. CT-based image fusion has been used to merge a Landsat enhanced thematic mapper plus, PAN and MS images. Based on experimental results using indicators of bias, the CT-based method provides better visual and quantitative results for RS fusion. The indicators are

- SD which indicates the dispersion degree between the gray values and the gray mean values,
- correlation coefficients (CCs) which indicate the degree of correlation between two images,
- the relative average spectral error characterizes the average performance of a method in the considered spectral bands,
- relative global dimensional synthesis error (ERGAS) which measures the spectral distortion between the reference image and the fused one.

Shearlet and nonsubsampled shearlet transform (NSST) are also used in this regard. In Ref. ^{69}, the authors propose in the first step to register two RS images. Then the shearlet transform is applied. The number of directions is set to 6 and scales to 5. The low frequency coefficients are selected based on the average rule while the high frequency coefficients are chosen using the maximum absolute value rule as mentioned in Eq. (28):

In Ref. ^{70}, the authors have used NSST.^{32} They have involved the SR^{71} paradigm in pansharpening application.^{72} The idea is about fusion the intensity component (I) of the MS image (the hue and saturation components are not used) with the PAN image. First, the authors decompose the images using NSST into low frequencies ${LI,LPAN}$ and high frequencies ${HI(s,d),HPAN(s,d)}$, where $s$ corresponds to scale and $d$ to direction. For the low frequency coefficients, they propose to construct vectors using patches from both $LI$ and $LPAN$. Those vectors are represented sparsely using a learned dictionary with K-SVD^{73} as mentioned in Eq. (29).

Then the new fused vector patch $vFi$ is reconstructed using the dictionary $D$ and the fused sparse coefficients are obtained as follows (31):

Otherwise, the high frequency coefficients are fused according to large local energy rule. This means that the energy of each scale $s$ and direction $d$ in the shearlet transform is calculated, as mentioned in Eq. (32) and then fused:

The experimental results show that the proposed method conserves better spatial and spectral information and is able to enhance the fused image better than Refs. ^{74} and ^{75}. The render result is illustrated in Fig. 16. We can notice how well preserved the spatial details and color information are. Fused images have better visual accuracy compared to existing methods: AIHS,^{76} AWLP,^{75} SVT,^{74} and ATWT.^{77}

The contourlets were also used^{78} in the same context to merge SPOT and ALSAT-2A images. The authors advanced a conception of a new fusion scheme by combining the PCA and the NSCT in order to overcome the spectral distortion caused by the PCA. The aim of this work was to find a compromise between enhancing the spatial resolution and preserving the spectral information at the same time. After decomposing the MS using PCA, a histogram matching is applied to adapt the contrast between PAN and the resulting components. Then the obtained components and PAN image are decomposed using NSCT into approximation coefficients app and details coefficients det. The fusion rules applied to this approach are represented by both Eqs. (33) and (34), where authors extract the approximation coefficients of the fused image $reconstim$ from the PAN image and the detail coefficients from the MS image.

Then the inverse of NSCT is calculated. The fused image obtained due to the proposed method, represented edges, contours of roads and buildings, and any structure shapes on the ground, better than other method, namely intensity hue saturation (IHS), PCA-IHS, and high pass filter (HPF). Moreover, the proposed method conserves the spectral information as well. Compared to PCA-NSCT based method, the resulting fusion image was not blurred.

Another injection model in pansharpening application is proposed in Ref. ^{79}. Authors worked on QuickBird and IKONOS-2 imagery. The injection model is built on an adaptive cross gain, i.e., a ratio of local SD. Both images are decomposed using curvelets and then merged together by applying interband structure model. Compared to IHS, ATWT, and HPF, the proposed method exceeds them since it has the best rate according to the indicators: ERGAS, spectral angle mapper, and Q4. The resulted fused image is visually superior and succeeds in producing a tradeoff between different sensors.

The authors of Ref. ^{80} propose an approach for multisensor image fusion, based on beyond the wavelet transform domain (CT, bandlet, contourlet, and wedgelet). The approach consisted of the following steps: first, the authors decomposed the images into coefficients using beyond wavelet transform. Second, they selected from the two images the low frequency coefficients using maximum local energy (MLE) rule. It is calculated in a local $3\xd73$ sliding window as shown in Eq. (36):

The use of several transforms in an approach could increase the quality of fusion. Works like Zhong et al.^{81} and Zhanga et al.^{82} combined wavelet and curvelet together, to enhance wavelets’ abilities.

Neural networks have been also investigated in RS fusion. MGD were combined with this powerful tool in order to provide more directionalities. In fact, in Ref. ^{83}, a fusion method of SAR images based on pulse couple neural network (PCNN)^{84} is proposed. The highlight of this algorithm is to use the global feature of source images and calculate the gradient information of them after being decomposed into several directions. Then wavelet is applied to decompose the images further while the high frequency coefficients are selected to be the input of PCNN to get a fire map. The max rule is applied to get the final fused image. The framework of the proposed method is detailed in Fig. 17.

PCNN was also used in Ref. ^{87}. The idea is to use the contourlet hidden Markov tree (CHMT) model^{88} to describe the statistical characteristics of contourlet coefficients of RS images. Besides, PCNN is used in this work in order to select the high-frequency directional sub-bands. First, contourlet transform is performed on the registered multisource SAR images and then the contourlet coefficients are trained by expectation maximization algorithm to calculate their edge PDFs. The highest magnitude in low-frequency sub-bands is selected. The high-frequency directional coefficients are updated by multiplying it by its edge PDF. Finally, a new clarity saliency measure is defined and used to fuse the high-frequency sub-bands. This method is compared to others: wavelets hidden Markov trees + PCNN and CT + PCNN. The results show that CHMT-PCNN can capture more directional information. Moreover, it needs less training time even with images abundant in textures. The criteria of assessment used in this work were the mutual information, weighted fusion (evaluates the visual effect), edge-dependent fusion quality index proposed by Ref. ^{89} (describes edge information in the image), and common information (represents the gradient information) proposed in Ref. ^{90}.