Unsupervised fuzzy-wavelet framework for coastal polynya detection in synthetic aperture radar images

The automated detection of coasts, riverbanks, and polynyas from synthetic aperture radar images is a difficult image processing task due to speckle noise. In this work we present a novel Fuzzy-Wavelet framework for bordeline region detection in SAR images. Our technique is based on a combination of Wavelet denoising and Fuzzy Logic which boost decision-making on noisy and poorly defined environments. Unlike most recent filtering-detection algorithms, we do not apply hypothesis tests (Wilcoxon-Mann Whitney-G0) to label the edge point candidates one by one, instead we construct a fuzzy map from wavelet denoised image and extract their borderline. We compare our algorithm performance with the popular Frost–Sobel approach and a version of Canny’s algorithm with data-dependent parameters, over a database of real polynyas and coastline simulated images under the multiplicative model. The experimental results are evaluated by comparing Pratt’s Figure of Merit index of edge map quality. In almost all test images our algorithm outperforms the standard algorithms in quality and speed. Subjects: Artificial Intelligence; Computer Science (General); Image Processing

ABOUT THE AUTHORS Karim Alejandra Nemer (S'01) was born in La Rioja, Argentina, on February 22, 1972. She received the degree in System Information Engineer from the National Technological University of Córdoba, Argentina, in 1996. In 2011 she became a PhD student in the National Technological University in Córdoba, under the advising of Prof. Ana Georgina Flesia (PhD), working on edge detection techniques for the unsupervised interpretation of synthetic aperture radar images using Wavelets, Fuzzy Logic, and Artificial Intelligence. She is interested in various applications of image processing using artificial intelligence methods, mathematical statistics estimation, probabilistic algorithms, and analysis of random processes. She is involved in several image processing applications related to industry, ecology, analysis of animal hair (camelids), medicine, and other fields. She works in collaboration with researchers from the FAMAF Universidad Nacional de Córdoba where PhD courses on image processing are given.

PUBLIC INTEREST STATEMENT
The study of the evolution of coastal areas with ice (polynyas) over time is essential to monitor climate and ecological changes and to analyze their repercussions on the planet. To make this study, satellites send microwave beams to these areas of the Earth and capture the reflections of the waves as synthetic aperture radar (SAR) images. These types of images are characterized by a noise, called speckled, which causes many miscalculations when the optical image processing is applied and leads to erroneous interpretations. Recently, several techniques that combine artificial intelligence mechanisms with signal analysis were implemented to overcome such drawbacks. This paper presents an image processing framework that combines the Wavelets Transforms with Fuzzy Logic to produce an homogeneous map of the image with few colors where the coastal polynyas edges can be trivially detected. This edge detection method is compared with two traditional edge detection methods giving excellent results.

Introduction
The study of polynya coasts is of paramount importance since it allows to analyze the effect of temperature variation and krill insertion in the ocean, among others (Stringer & Groves, 1991;Tamura & Ohshima, 2011). An accurate determination of the coastline contributes to monitor its evolution over time. The images taken of these zones are mainly of synthetic aperture radar (SAR) type due to climatic reasons. Satellite born SAR images represent the complex reflectivity map of a scene, and unlike optical images, the interpretation is not consistent with common visual perception. Furthermore, the direct application of conventional image processing tools conceived from an optical point of view usually gives suboptimum results on SAR data (Flesia, Gimenez, & Rufeil Fiori, 2013;Lim & Jang, 2002). A crucial point for SAR image automatic interpretation is the low level step of scene segmentation, i.e. the decomposition of the image in a tessellation of uniform areas. Edgebased segmentation schemes aim to find out the transitions between uniform areas, rather than directly identifying those areas. The related algorithms generally work in two stages: they firstly compute an edge strength map of the scene and finally extract the local maxima of this map. The first stage can be achieved by filtering, creating well defined homogeneous areas (Alonso, López-Martínez, Mallorquí, & Salembier, 2011), or by using an statistical edge detection filter, which selects few points in the boundary of the object or section to underline (Girón, Frery, & Cribari-Neto, 2012). Popular filters that made use of statistical properties of SAR data are Frost and Lee filters (Frost, Stiles, Shanmugan, & Holtzman, 1982;Lee, Wen, Ainsworth, Chen, & Chen, 2009); more recent approaches apply wavelet shrinkage to generate a skeleton edge map. Testing methods consist in scanning the image with a sliding two-region window and evaluating for each position if there is a change between the two regions by hypothesis testing. Such tests are dependent of the true underlying conditional distributions that are emitted by the states (classes) of the image field, and the regions considered in the test implementation (Baselice & Ferraioli, 2012;Lim & Jang, 2002;Touzi, Lopes, & Bousquet, 1988). In all cases a standard edge detector like Sobel, Canny, or Snakes, is applied to the skeleton intermediate image to obtain a map with single pixel curves (Gambini, Mejail, Jacobo-Berlles, & Frery, 2006). Alonso et al. (2011) advocate for the use of the bidimensional undecimated discrete wavelet transform (2dUDWT) for filtering noise while enhancing edges on SAR images. They propose a deterministic multiscale thresholding scheme on the Haar wavelet transform to produce an edge strength image that single out edge candidates. The final map is made by applying a snake algorithm to the edge strength image to get one pixel wide continuous edges. These combination algorithms are popular in the literature of SAR edge detection. Radiometric edge strength masks are usually obtained by hypothesis testing, using the likelihood ratio filter, the robust rank order or the Wilcoxon test (Dias, Cribari-Neto, & Ospina, 2015;Girón et al., 2012), or costume made maximum likelihood ratio tests (Gambini, Mejail, Jacobo-Berlles, & Frery, 2008;Gambini et al., 2006;Lim, 2006). The main drawbacks of such radiometric approaches lie on the slow processing speed, and the need of hypothesis about the distribution of the SAR region where the borderline is sought. Non parametric approaches bypass the need of knowledge about such distribution but require more data on the sample. Moreover, spatial information is lost, since the pixel samples are considered independent and identically distributed. Flesia et al. (2013) discussed the idea of applying a Hidden Markov Model to the coefficients of the bidimensional Haar UDWT, selecting edge pixels if their strength was propagated over several bands. Such approach reduce the distribution of information requirement of former methods, but still considers the coefficients as an independent identically distributed sample that needs to be merged to generate the edge strength map, losing all spatial information.
To overcome these constraints, we propose a new algorithm that projects and reconstructs the image using a one dimensional Discrete Wavelet Transform (1d-DWT) on rows and columns of the image independently, where the reconstruction uses only information from the J-band of the 1d-DWT. After reconstruction, both images are merged using Fuzzy Logic. The result is a smooth map with sections related to the sets considered in the Fuzzy Logic. The borderline of such sections is extracted by applying the Canny's algorithm, generating the final edge map.
Our approach is non-parametric and faster than the ones introduced by Flesia et al. (2013) and Alonso et al. (2011), which rely on the Haar UDWT. Besides that, our algorithm is flexible in the sense that it allows the use of different mother wavelets for wavelet analysis and reconstruction, when those proposals are based only on Haar wavelet analysis. In our experiments, we tested our algorithm with the Daubechies (2, 3, 4, 5 and 6), Haar, Symlets, Coiflets, Discrete Meyer, and Reverse Biorthogonal wavelet families, being the Discrete Meyer wavelet which yield the most accurate and robust results. Moreover, Flesia et al. (2013) and Alonso et al. (2011) merge the wavelet related information using simple multiplication (logic AND) and addition (logic OR), which are hard deterministic rules. Using Fuzzy Logic it is possible to obtain continuous levels of membership (between 0 and 1) to a given set, since it is based on heuristics rules of the type If (antecedent), then (consequent) (Driankov, Hellendoorn, & Reinfrank, 1996). Within the Artificial Intelligence, Fuzzy Logic is used to solve a variety of problems, and involves decision-making using incomplete, vague, and even with some level of error, data. Fuzzy logic has been applied on previous work on image segmentation and classification of polynyas on SAR imagery (Devi & Veera Senthil, 2013;Klir & Yuan, 1995). In this paper we aim to develop an automated edge and borderline detection, but our intermediate product, which is a Fuzzy Logic segmentation, may be shaped into a fully segmentation outcome. We will study its possibilities elsewhere.
Among the most important differences between optical and SAR images are the appropriate noise models (Girón et al., 2012). Moreover a complex SAR image may be represented as the convolution of the local complex reflectivity of the observed area with the impulse response of the SAR system. Hence, the SAR data pixels are roughly the low pass filtered version of the complex local scattering properties of the observed scene. A usual model for the complex image is obtained considering its real and imaginary parts as independent, zero mean and equally distributed random variables. Consequently, the intensity of a SAR image follows a negative exponential distribution, and can be further modeled as a product of the radar cross section and a noise term called speckle, which is exponentially distributed with mean equal to one (Alonso et al., 2011). Data format (complex, intensity or amplitude) and polarization are important information to consider in simulation, since the statistical model of data must be changed accordingly. In this paper we design a careful simulation of SAR coastline data under the multiplicative model of SAR monopolarized intensity data, following Frery, Muller, Yanasse, and Sant'Anna (1997) and Flesia, Lucini, and Perez (2015), and the robust edge detection framework proposed by Lim (2006). This paper is organized as follows: In Section 2 we present the mathematical background of the DWT and Fuzzy Logic. The details of our method is explained in Section 3. In Section 4. We evaluate the results of our method using two filter + edge standard detector approaches (Frost-Sobel) and Canny's method with data-dependent estimated parameters, over a data-set of synthetic images and real SAR images of polynyas (ASF, 2014). The synthetic example was made by simulation of the G 0 SAR image model, which divides image regions on three sets according to texture: extremely heterogeneous, heterogeneous and homogeneous (Flesia et al., 2015;Girón et al., 2012). Accuracy was measured with Pratt's Figure of Merit (PFoM), using the stencil as ground truth in the simulated example and handmade ground truths in the polynya database; an effective method must produce an average PFoM greater than 0.8 (Pratt, 2007). The last section was left for conclusions and prospects.

General background
In this section we describe the theoretical background necessary to understand the detection problem and to describe our novel method. This includes Polynyas, Wavelet Transforms, and the Fuzzy Logic decision-making process.

Polynyas
The accurate detection of coastal polynyas is essential for the estimation of sea ice production (Tamura & Ohshima, 2011). When the polynya is limited on one side by the coast, it is called a shore polynya, and when it is limited by land-fast ice it is called a flaw polynya. If it recurs in the same position every year, it is called a recurring polynya. Figure 1 shows SAR images of polynya coasts of different zones of the Antarctic. In some parts of the images, the limit between the water (black) and the ice (clear area) is easily seen on an equalized image, but in some other parts, the ice is so thin that the limit is difficult to determine. SAR imagery is not made for optical discrimination, thus statistical discrimination between regions is customary. The multiplicative model (Oliver & Quegan, 2004) is a description of the data that has enabled the best statistical analysis discussed in the literature. Within this model, texture is characterized by statistical distributions related to the G 0 distribution (Frery et al., 1997), separated in three categories, homogeneous, heterogeneous, and extremely heterogeneous. Speckle noise contaminates the image with a large variance distribution, but the texture is given by the image, the heterogeneity and granularity of the regions are not due by noise, but they are truly related to the scattering process and characterize the regions.

Discrete Wavelet Transform (DWT)
Wavelets are localized basis functions, so they are good for representing short-time events. More precisely a wavelet is a short-duration oscillating signal whose energy is finite and is concentrated in a given time interval. They satisfy specific mathematical requirements and are used to represent data or other functions (Mallat, 2008). A Wavelet Transform WT may be loosely defined as a projection into a wavelet basis or frame that generates information blocks (coefficients) in scaling and time from one signal, by translation and dilation operations of a single fixed function called mother wavelet ψ(t) defined as in Equation (1): where a and b are real numbers, a allows the signal expansions and contractions, and b allows the change of the signal position over time (in the context of this work, instead of time t, the wavelet is evaluated in spatial form using x and y as variables on the rows and columns of image pixels). Wavelet transforms typically involve waveforms that have finite support, that is, are nonzero on a finite length interval, in contrast to cosines and sines that oscillate for all time. Whereas the discrete Fourier transform (DFT) decomposes a signal x into a sum of sinusoids at particular frequencies, the DWT splits x into components called wavelet details and smooths that are associated with different frequency bands.
In signal processing, a filter bank is an array of band-pass filters that separates the input signal into multiple components, each one carrying a single frequency sub-band of the original signal. Mallat's pyramid algorithm for multiresolution connects wavelets to filter banks; the WT coefficients at each scale are filtered and sub-sampled to give coefficients at the next scale, starting with the sampled signal. This is accomplished by repeatedly applying a High-Pass wavelet Filter to generate wavelet coefficients W j, k and a Low-Pass scaling Filter to generate scaling coefficients V j, k . See Figure 2 for graphs of the wavelet and scaling filters used in this article.
(1) The wavelet coefficients W j, k yield the level j wavelet detail D j , and the scaling coefficients V j, k yield the level j wavelet smooth signal S j , via application of the inverse DWT. The original signal x is the sum of the wavelet details D 1 , …, D J plus the final wavelet smooth S J , where J is the number of levels computed. Figure 2 illustrates this process for J = 3. The wavelet detail D j corresponds to a frequency band with period range 2 j Δt to 2 j+1 Δt time units, where the signal x was created by sampling every Δt time units. Hence, each level is associated with a scale twice as long as the previous level s. The wavelet smooth S J is obtained by repeated smoothing via the LPF and so forms the trend of the signal. Together, the wavelet details and smooth form a multiresolution analysis of the signal. Since a downsample operator is applied after each filter, the number of coefficients match the number of samples of the signal.

Fuzzy Logic decision-making process
Fuzzy Logic is a widely used process, through which it is possible to obtain continuous levels of membership (between 0 and 1) to a given set. It is based on heuristics rules of the type If (antecedent), then (consequent) (Driankov et al., 1996). Within the Artificial Intelligence, Fuzzy Logic is used to solve a variety of problems, and involves decision-making using incomplete, vague, and even with some level of error, data. Figure 3 shows the membership of 12 elements to three fuzzy sets. It can be seen that there are elements belonging to two of the sets; it is also possible that elements belong to all the fuzzy sets, with certain levels of membership. With conventional logic, computers can manipulate strictly dichotomic values, such as true/false or yes/No. In Fuzzy Logic, mathematical models are used for representing subjective notions such as hot/warm/cold, that are converted into real values between (0,1), so they can be manipulated by computers (Senthilkumaran & Rajesh, 2009). This means that the quality quantifiers are translated, within the language, and that they are transformed in continuous numeric values.
Mathematically, we say that the element x belongs to A, x ∈ A, with a number μ A (x) ∈ [0, 1] that stands for the degree of membership of element (x) to the set A. If μ A (x) = 0, then x ∉ A. If μ A (x) = 1, then x belongs entirely to A. If 0 < μ A (x) < 1, then x belongs partially to A. This stands for elements as well as for sets, this is why the B set may belong partially to the set A. Fuzzy Set Theory includes traditional operations, such as union, intersection, difference, negation and complement. Also, there are concepts belonging exclusively to fuzzy sets: • The nucleus of a fuzzy set is comprised of all elements belonging entirely to said set, i.e. μ A (x) = 1.
• The support of a fuzzy set is the set of elements belonging to some extent the fuzzy set, i.e. μ A (x) ≥ 0.
• Inclusion: Let A and B be two fuzzy parts of the C set. It is said that i.e. all elements of A belong mostly to B. In Artificial Intelligence, Fuzzy Logic is used for decision-making, text comprehension and for understanding of information. The systems that use it are generally more robust, and tolerant to noise and lack of data (Klir & Yuan, 1995).

The proposed method: FWF
The borderline detection technique to be developed in this section is based on the theory described in the previous section. The approach that was followed was to first derive a smooth collection of vectors to estimate the strongest discontinuities of the textured backscatter X(x, y) from the observed image Z(x, y) under the assumptions of stationary image data. The non-stationary aspects of real SAR data are then incorporated by identifying the degree of membership of such discontinuities that vary with position. The result of this analysis is an adaptive technique which is computationally efficient, provides a texture segmentation estimate in coherent areas, and tends to preserve the edge structure of the image.
The fuzzy wavelet framework (FWF) consists in: Step 1: Denoising by wavelet thresholding (Figure 4). For each pixel we apply the Discrete Wavelet transform (DWT) with preset level of decomposition (J) in rows and columns, independently, to increase speed and precision. After this, we apply the inverse wavelet transform, keeping only the last level of decomposition, setting all other coefficients zero. Here two intermediate images are obtained, one for rows and other for column processing. The objective of this step is to remove the noise, and obtain a lower level noise image, to be processed in the following stage. To perform this step the level of decomposition J must be selected.
Step 2: Merging step ( Figure 5 illustrates this and next steps). We combine linearly the pixels of both the intermediates images, of all the pixels in the fuzzy area. To perform this step a neighborhood of pixels must be selected and defined as the fuzzy area radius (FAR), being a neighborhood of 1 for FAR = 1, of 9 for FAR = 2, and of 25 for FAR = 3.
Step 3: Membership generation. We apply a nonlinear operator to the merged map. In this case, the operator is the Probabilistic-Or (P-Or), as presented by Devi and Veera Senthil (2013), obtaining a fuzzy map of outputs in which the level of membership of each pixel to each segment is shown.
Step 4: Defuzzification step. For image defuzzification, a Bayesian criterion is applied, taking the greater probability of membership of each pixel. In the defuzzification process, the candidate pixel is taken with a determined neighborhood defined by the FAR coefficient. The number of sections in the image must be specified, as this determines the fuzzy sets onto which the operator is applied. Here a final coherent map with S-sections is obtained.
Step 5: Borderline detection. We identify borderlines as boundaries of the sections map, generating the final edge map that will have the value 0 for non edge pixel and the value 1 for an edge pixel.

Characteristics of the wavelet used
The DWT is a flexible tool, since allows the selection of the appropriate type of wavelet for detecting specific features within the image. In our simulated experiments we studied 10 types of wavelets (modifying slightly the algorithm accordingly when changing from orthogonal wavelets to other types of wavelets). We obtained the best results with the discretization of the Meyer wavelet. This wavelet and its scaling function are defined in the frequency domain; the wavelet has not compact support, due to his good localization in frequency, but decays to zero faster than any polynomial, and it ensures orthogonal analysis. Although the Meyer wavelet is not compactly supported, there exists a good approximation leading to FIR (Finite Impulse Response) filters which can be used in the DWT, called Discrete Meyer. With this particular wavelet we obtained a mean 15% increase in PFoM accuracy on the whole simulated experiment.
Wavelet analysis has a fundamental factor in the Level of Decomposition (J) which has to determined in advance. Since we are also thresholding by keeping only such level of coefficients, the decision become even more important in the images with extremely heterogeneous texture. The DWT is a dyadic transform; for instance, J = 3 implies the use of eight point samples to produce a coefficient. The Frost filter that is used in the comparison section is adaptive, the number of samples included in the smoothness step is determined by the variance of the image. Figure 6 shows a simulated image with two extremely heterogeneous texture (a) and its ground truth (b), a coastline one pixel wide simulated with a random walk. The output edge map for each parameter J = 1, 2, 3 can be observed in Figure 6(c), (d), and (e).
In this figure we observe the tradeoff between high denoising and detail. The larger the level of decomposition, the fewer noise is taken as edge point, although the resulting borderline is smoother, losing part of the detail. Denoising thus should not be overdone, since fuzzy membership will also remove texture related pixels that could remain in the intermediate map.

Options for FAR
The merging and membership labeling steps of the method need the determination of the FAR and the merging nonlinear operator. In this work, the function used to determine membership, i.e. the set to which each pixel belongs, is the Probabilistic Or (P-Or). Such nonlinear operator responds to the function given by the Equation (2): where 0 ≤ a, b ≤ 1, and graphically behaves as shown in Figure 7. This operation is performed over both intermediate images generated by the Wavelet stage, that is, on intermediate columns and rows images.
The neighborhood on which the P-Or function is applied is another key factor for the proper segmentation. When the image texture is homogeneous, only the evaluated pixel is used (fuzzy area radius FAR = 1), while in the cases of heterogeneous texture, a neighborhood of 9 pixels is required (fuzzy area radius FAR = 2), and in the case of extremely heterogeneous texture, a neighborhood of 25 pixels is needed (fuzzy area radius FAR = 3).  (e) show the edge outputs after applying the method with a FAR equal to 1, 2, and 3. As with the decomposition level, the larger the radius of the diffuse area, the fewer texture pixels are taken as edge pixels and only the region's borderline remains. Nevertheless, the borderline is smoother than the Ground Truth, losing part of the detail.

Results and discussion
In this section our primary concern is to compare our detector with two existing filtering-detecting proposals such as the classical Frost-Sobel detector (Dias et al., 2015;Frost et al., 1982) and the Canny detector (Canny, 1986) with data-dependent parameters. For such task, simulated and real images processed by the three detectors are presented and evaluated in this section. This assessment is performed in two ways. First, a quantitative evaluation of the algorithms is made with the edge Figure of Merit proposed by Pratt (2007). Then a qualitative evaluation is conducted using real sensor images; in this case, data of coastlines and polynyas adquired by the Sentinel-1 satellite.

PFoM evaluation of performance on synthetic images
The filtering/detecting algorithms are evaluated here using the following procedure (Frost et al., 1982): • produce the simulated input images, • estimate the parameters of filtering and threshold operators, • apply the enhancement operators (wavelet denoising, adaptive filtering, Gaussian filtering), Figure 5. Merge: Reconstruction process. Using both intermediates images, a smoothing filter with the FAR coefficient is applied. After that, we apply a non linear operator and then a defuzzification process to calculate the section coefficient for each pixel, to obtain the final homogeneous map. Finally, we apply an edge identifier to obtain the final edge map.
• apply Sobel, Prewitt, or Robert or any gradient operator to each enhanced image, • apply a threshold to the gradient image, and • apply PFoM index.
We are comparing our method with two filtering-detecting algorithms, Frost-Sobel and Canny (Canny, 1986;Frost et al., 1982). Frost method estimates the parameters in an automated fashion from data. Canny method applies a Gaussian filter with parameter σ, a data-dependent thresholding, T 1 and T 2 , of the gradient image and a hysteresis step that links isolated edges, producing single pixel edges. We conducted a simulation experiment to determine the best parameters for Canny's method. We decided that a data-dependent formulae was needed to obtain good results, thus the σ coefficient is computed for the different types of texture, as shown in Table 1, and depends on the standard deviation S D of the images, the number of sections S in the images, and M, the mean of the grey values of the pixels in the image. Using these parameters, the low and high thresholds can be computed as in Equation (3): Parameter determination is a subject of research itself (Gimenez, Martinez, & Flesia, 2014). Canny's parameters are not automatic, in the sense that they depend on the number of classes in the image and several constants that have been selected by trial and error, searching for the best performance over the whole database, which depend of the level of texture of the images. For the Sobel method, we use Matlab automatic cutoff based on a root mean square estimate of noise; the mean of the magnitude squared image is a value that is roughly proportional to the signal to noise ratio. This cutoff value does not take into account the number of sections (regions) in the image or the degree of texture. We set this parameter in this way on purpose, since usually there is no knowledge about the number of classes in the image, or the gray levels standard deviation given the texture class. Results from Subsects 4.1.3 to 4.1.6 will show the performance of these filtering/detecting algorithms implemented with prior knowledge and without it.

PFoM
In order to evaluate the performance of the detection, we use the PFoM. This index considers the following information to score the accuracy on edge points. Pratt (2007) states that there are three types of errors associated with edge detection: (1) Missing valid edge points, (2) Failure to localize edge points, (3) Classification of noise fluctuation as edge points.
The PFoM index balances the three types of error, and it is defined by Equation (4): (3) • Get the rows and columns of each Ground Truth edge pixel and save them in a matrix of two columns, the first for the row and the second for the column of each edge pixel. The number of rows of this matrix is the number of edge pixels in the GT.
• For each edge pixel of the actual edge map image of the wavelet stage, get the Euclidean distance to all the edge pixels of the GT, take the minimum value and save it in a vector of distances d. The d vector has as many elements as number of pixels are in the evaluated edge map image.
• Compute Equation (5) using the values stored in d.
The scaling factor used is a = 1 / 9, as used by (Pratt, 2007). PFoM index images with values from 0 to 1, being 1 the optimum response. If the response is close to 0, the map has the lowest quality.

Database and FWF parameters
In what follows, we will report the results of several Monte Carlo simulations, which were performed to assess the finite sample merits of all borderline estimates in SAR images. All simulations were carried out using Matlab programming language. The hardware used was an Intel Core i7 CPU 2.4 MHz 8G RAM computer running on Ubuntu 12. The simulated data have been generated according to the statistical distributions stated in Table 2, which correspond to the multiplicative model of SAR imagery proposed by Frery et al. (1997). According to this model, the data are described by a random variable Z which can be viewed as the product of independent random variables X and Y, Figure 7. Probabilistic-Or operator applied to two datasets sampled from different statistical distributions. When working with imagery that has texture strongly related with high variance statistical distributions, Fuzzy Logic helps the determination of segmentation thresholds. where X models the properties of the imaged area (backscatter) and Y models the multiplicative noise (speckle) as shown in Equation (5): where (s 1 , s 2 ) denotes the spatial position of the pixel. This model has been widely used in the processing and analysis of SAR images (Dias et al., 2015;Flesia et al., 2015;Frery et al., 1997).
A set of 400 basic template images of 2 and 3 and sectors (200 of each) were generated with a ground truth following a random walk. Most edge borderline simulation experiments are made with straight edges separating two sections, as used by Dias et al. (2015), which is a simplistic situation good for studying performance of algorithms based on two-population statistical tests. In this paper we aim for a more realistic situation with respect to the GT generated (a random walk output).
Therefore, a database of 4,800 specked images were generated using 12 types of images characterized by the texture statistical distribution for 2 and for 3 sections. To perform the edge detection, 10 different Wavelets with 1 to 5 levels of decomposition were applied to each of these noisy images. These results were processed with the P-Or function with FAR of 1, 2, and 3. A total of 1,440,000 results were obtained.
A preliminar analysis of these results leads to conclude that the factors which have more influence on the PFoM index are: • Number of sections in the image, S.
• Level of decomposition, J.
The rest of the factors have less influence on the performance of the FWF. In particular, by changing the type of wavelet, the PFoM results show a mean increase of 15% when using Discrete Meyer wavelets. Given the lack of space, we will report only the outputs of our method with it, so, when we talk about FWF, we refer to the Discrete Meyer Wavelet exclusively.
In the following subsections, we will discuss in detail the incidence in the PFoM values of the aforementioned factors, and we will report the set of parameters, pair (J, FAR), that gives better performance. Global performance will be assessed in two different ways: reporting relative improvement using Canny and Frost as baseline methods, and p-values of nonparametric Friedman's test to compare methods effects in a two way layout.

Homogeneous intensity backscatter images
The distribution of speckled homogeneous textures under the multiplicative model (Frery et al., 1997) follows the Equation (6): Z I (s 1 , s 2 ) ∼ Γ(n, n∕ 2 ) The Monte Carlo experiment for this type of texture set de step of the edge Δn to 40 or 60 and Δ as 4 or 8. The number of Monte Carlo replications is R = 200, for each of the four templets. After evaluating the performance of our method with parameters J ranging from 1 to 5 and FAR ranging from 1 to 3, we report that (J, FAR) = (1, 1) is enough to provide an average PFoM of 0.9548, with a standard deviation of 0.0015. Figure 9(a) shows an example of a synthetic image generated with Equation (5); Figure 9(b) shows the GT of this image. Figure 9(c) shows the result of applying Canny's, Figure 9(d) shows the result of applying Frost-Sobel and Figure 9(e) shows the result of applying FWF. We thus recommend this set of parameters for borderline discrimination of regions with homogeneous texture. We should notice that the modes of the distributions of the sections are fairly separated, thus all methods should have a very good performance despite the high variance. This result is used as baseline test to check if the algorithm implementation is correct.

Heterogeneous intensity backscatter images
The density that characterizes this type of texture is given in Equation (7) where n = (40, 60), α = (4, 6) and = (3, 5), (Frery et al., 1997). After evaluating the performance of the FWF method with parameters J ranging from 1 to 5 and FAR ranging from 1 to 3, we report that (J, FAR) = (2, 3) has the best performance. Edge detection in this example is more challenging than the previous one because the edges divide the image in regions with textures that are visually quite diffused. The detection results are displayed in Figure 10. It shows an example of an artificial image with noise and backscatter characterized by the indicated coefficients, and the result of applying the FWT with Discrete Meyer Wavelet, J = 2, FAR = 3. A PFoM = 0.9123 was obtained for this example.

Extremely-heterogeneous intensity backscatter images
The function that characterizes this type of texture is (8) where n = (40, 60), = ( −1 5, −8), γ = (1, 8), (Frery et al., 1997). After evaluating the performance of our method with parameters J ranging from 1 to 5 and FAR ranging from 1 to 3, we report that for this type of image, three levels of decomposition and a FAR of 3 have the best performance. Figure  11 shows an example of artificial image with noise and backscatter characterized by the indicated coefficients, and the result applying the FWT with Meyer Wavelet, Level of decomposition J = 3, FAR = 3. A PFoM = 0.9013 was obtained.

Global performance
The performance of FWF, Frost-Sobel, and Canny's method, for synthetic images corresponding to all types of backscatter and for both 2 and 3 sections analyzed can be seen summarized in Figure 12.
The FWF parameters were set as reported before for each type of texture. It can be seen that the FWF has better global performance than the other methods, but also that the values of PFoM are  , n) high, FWF has PFoM over 0.8 on more than 90% of the cases and 50% of the cases is better than 0. 94. Frost/Sobel with automatic parameters is the worse of the methods, since it has not incorporated the idea of borderline edge. Frost filter is too gentle, leaving texture in the image and Sobel picks small runs of texture points as small edges. These mistakes are heavily penalized by PFoM. With Canny we have chosen a smoothness parameter high enough to smooth out the image with a large variance filter when needed, and that extra knowledge was of paramount importance.
The distributions of the simulated imagery under the multiplicative model have large variance, therefore the image factor introduce a large source of variation in the possible values of PFoM. Two reduce the incidence of such factor, and to report the degree of improvement of FWF relative to the baseline method, we introduce a new index in terms of the PFoM as follows, where GT indicates the Ground Truth image, I is an edge map, and sub-indexes a and b are the tested and basis methods through which the edge map images I a and I b were respectively obtained. To define this relative index we found inspiration on the relative improvement index proposed by Joossens and Croux (2004) and Fraiman, Justel, and Svarc (2008) for reporting improvement on classification and clustering rates.
Using this new index, we score the complete image data-set for 2 and 3 sections to compare FWF with Canny as baseline method. A histogram of such values is shown in Figure 13(a) and (b) for 2 and 3 sections respectively. We also report the values of the relative index for FWF with Frost-Sobel as baseline method in Figures 13(c) and (d) for 2 and 3 sections respectively.
Each histogram shows the frequency of the dPFoM values obtained with the FWF relative to the baseline method. When the dPFoM is positive, the result is better for FWF, otherwise it is worse than the baseline method.
The Figure 13(a) and (b) shows that the FWF obtains better PFoMs than Canny, in almost all the studied cases, except for a few cases in heterogeneous backscatter. In average, the relative improvement of FWF over Canny is 5.60% for 2 sections, and 4.62% for 3 sections. Relative improvement climbs to 10% on heterogeneous textures on both two and three sections.
The Figures 13(c) and (d) shows that the FWF obtains better PFoMs than Frost-Sobel, in all the studied cases. In average, the relative improvement of the FWF over Frost-Sobel is 57.13% for 2 sections and 38.41% for 3 sections.
To determine whether the results between different methods have statistically significant differences we have performed the Friedman's Analysis Of Variance (ANOVA) test, for all three methods, eliminating the image factor. The Friedman test is a nonparametric test that compares three or more matched or paired groups, which is used when the same parameter has been measured under different conditions on the same subjects. In our case this test evaluates the methods and eliminates the image factor, which is very relevant in a nonparametric Friedman's study. In the result of the test, one can see that the null hypothesis of equal mean ranks is rejected with a p-value smaller than the precision of Matlab. Therefore, there is a significant difference in quality in the edge maps. The relative improvement histograms of Figure 13 depict the extent of the difference betweeen the methods, measured by PFoM. Figure 12 show the high PFoM scores of FWF.
Friedman's test indicates that there is significant difference between all the methods, but this is not sufficient to determine whether the proposed method is significantly different to each of the other two, Canny and Frost-Sobel. To perform this analysis, the results provided by the Friedman's test were taken to conduct a multiple comparison of mean rank with all the three methods, shown in   Figure 14(a). The same sequence of two analyses was performed for FWF and Canny methods, as shown in Figure 15(b), for FWF and Frost-Sobel, as shown in Figure 15(c), and for Cannny and Frost-Sobel, as is depicted in Figure 14(d).
We conclude from this findings that FWF is significantly different from Canny and Frost-Sobel, either taken all together or by pairs.

Real SAR images of coastal polynyas
SAR images of polynyas were used to test the method, obtaining a level of PFoM ≈ 0.8. Part of this error is attributed to the GT miscalculation. This is shown in the images in the first column of the We should notice that FWF has PFoM over 0.8 on more than 90% of the cases and 50% of the cases is better than 0.94. Canny's method has more than 75% of the cases a PFoM larger than 0.8, but the best performance is 0.92. We will see in the next section than Frost-Sobel is good only with homogeneous texture. Figure 15. It is difficult to determine if the accuracy of the GT is of 100%, which allows us to infer that a PFoM of 0.8 is satisfactory.
It can be seen that the proposed method segments with good precision the coastlines of the polynyas. This enables us to label properly the kind of soil present in the image.

Conclusions
Many environmental studies involve the accurate determination of borderlines, like coastlines and river banks, from imagery taken with different sensors, atmospheric conditions and resolution through time. In the case of SAR imagery, it is also subject to speckle noise which is of multiplicative nature, increasing drastically the variance of the data when the region observed is considered extremely-heterogeneous. In particular, the study of coastal polynyas is involved in many environmental applications such as predicting and understanding climate change, and tracking the highest concentrations of phytoplankton, the foundation of the marine food chain.
The aim of this paper is to introduce a new method for the extraction of borderlines based on Fuzzy Logic and Wavelets that can accurately extract thin (one pixel wide) continuous curves in the position of change of the texture statistics of a SAR image defined under the multiplicative model. The method has a denoising step performed by wavelet hard shrinkage and merging and cleaning steps made using Fuzzy Logic. The degree of smoothness of the curve can be selected by the level of the wavelet decomposition and FAR. The Fuzzy Logic nonlinear merging operator is a Probabilistic OR, which helps cleaning the remains of true texture in the denoised image. In the SAR edge detection literature, there are two main approaches and both have statistics of the image into account. The first is a filtering-detection scheme and the second is a edge strength detection by hypothesis testing followed by detection. Our method belongs to the first category. Most algorithms of the second category are slow and do not perform well on jagged contours, as polynyas and other coastlines.
In this paper we found that when detecting borderline between simulated texture under the SAR multiplicative model, wavelet Fuzzy Logic approach provide substantially better detectors than Frost-Sobel and an enhanced Canny, for all the types of texture. The accuracy of the detection is