An Adaptation of Fourier Descriptors’ Spectral Energy for the Classification of Skin Cancer Melanoma

Skin cancers are the most common cancer throughout the world. In 2022, approximately 1.5 million new cases were reported. Approximately 330,000 of these new cases are melanoma. The Global Cancer Observatory estimates 331,722 new cases of melanoma and 58,667 deaths worldwide. Globally, melanoma comprises only 20% of all skin cancers, but causes the majority of skin cancer deaths. The corresponding mortality rates are about 0.53 per 100,000, respectively. Incidence is much higher in fair-skinned populations. For example, Australia, New Zealand, and Western Europe have the highest rates compared to Asia and Africa. An estimated 1.26 million people have been diagnosed with melanoma in the past five years, predominantly in Europe and North America. Melanoma occurs more frequently in men than in women in most regions. As an illustration of the geographic burden, Europe and North America together account for 78% of melanoma cases and 60% of deaths [1].

Skin cancer also imposes a heavy economic burden. Detailed global cost data are scarce, but country-level analyzes highlight the scale. In the United States, for example, nearly 6 million people are treated annually for skin cancer, at an annual medical cost of about 9 billion USD. Melanoma in particular causes substantial loss of productivity from premature deaths and morbidity. A US review estimated that melanoma mortality alone (and its associated years of life lost) corresponds to on the order of 2 to 3 billion dollars per year in lost productivity. Extrapolating to the global scale (with many more cases and much higher treatment costs in absolute terms) implies a global economic toll that is likely to be in the tens of billion annually. These figures do not count out-of-non-medicals or non-medical impacts of skin cancer [1].

Skin cancer begins when skin cells undergo changes in their DNA. This DNA contains instructions that guide how cells function. In normal cells, these instructions regulate growth, division, and when cells should naturally die. However, in cancer cells, DNA mutations alter these signals, causing cells to grow and divide uncontrollably. Unlike healthy cells, cancer cells don't die when they should, leading to an excessive buildup of cells. While the exact cause of skin cancer is not always known, it is generally linked to exposure to ultraviolet (UV) radiation from the sun. UV radiation can damage the DNA in skin cells, which may trigger cancer development. Globally, cutaneous malignant melanoma form of skin cancer occurs more frequently in people with fair skin.

Visual observation of a lesion is the most common method because melanoma lesions have non-symmetrical and irregular shapes, multicolor in different regions and larger sizes. However, examination of a physician with the naked eye may affect examination performance negatively because for some cases there are similarities, and differences may not be noticeable by the naked eye. It is reported that using computer-aided diagnostic tools improves diagnostic sensitivity by 20–30% compared with the examination of a doctor with the naked eye [2]. Therefore, many studies have been carried out using machine learning (ML) and deep learning (DL) techniques to diagnose the melanoma exceeding the performance of experienced dermatologists in detecting and diagnosing skin lesions. ML techniques include Decision Trees (DT) [3], Support Vector Machines (SVMs) [4], K-Nearest Neighbour (KNN) [5] and Artificial Neural Networks (ANNs) [6]. DL techniques include Convolutional Neural Networks (CNNs or ConvNet) and pre-trained models such as GoogLeNet, InceptionV3, ResNet, SqueezeNet, DarkNet, DenseNet, Xception, Inception-ResNet, Nasnet, and EfficientNet. Some of these studies have been examined by Grignaffini et al. with a systematic review [7].

In this paper, we contribute a computer-aided decision support system to identify melanoma skin cancer. We have extracted a novel feature from lesion image Fourier Descriptor’s (FD) energy spectrum using Ordinary Least Squares (OLS) method filtering out low frequency components of FDs using PCA + HMM. SHAP (Shapley Additive Explanations) analysis has proved the effectiveness of this feature [8]. We have included average color and shape parameters such as major/minor axis, eccentricity, compactness and perimeter as additional features. This way, we have minimized the number of features obtaining a moderate size ANN model. High accuracy with short train time and low processing power have been achieved compared to the models being used in the previous researches in literature. From the point of model complexity, training load and inference time, we observe that our ANN model has many advantages over recent deep learning models.

The organization of this paper is as follows: In the Related Works section, we will provide an overview of literature survey classified as machine learning and deep learning studies. In the Materials and Methodology section, we provide general information about how we selected dataset from the International Skin Imaging Collaboration (ISIC) database [9] as well as the cleaning of the lesion images in the Dataset Preparation subsection. In the Feature Extraction subsection, a novel approach will be introduced to extract a feature using Fourier Descriptors energy spectrum. ANN model details will be provided in this section as well. Our research’s performance results will be provided compared with previous studies in the Experimental Results and Analysis Section. In the Conclusion Section we will give the final suggestions as well as providing a guide for future research.

In this section, we will first introduce the dataset where we obtained our lesion images. Next, we will explain how the novel feature is extracted. After the feature extraction we will present our classification method based on machine learning model: ANN. Figure 1 summarizes the method used in this study.

Figure 1. Methodology of this study

3.1 Dataset preparation

Total of 1000 images used in this study were downloaded from ISIC (The International Skin Imaging Collaboration) website [9]. Classes were not balanced in ISIC’s datasets. Majority of the lesions were labeled benign. In the ISIC’s dataset, while most of the melanoma lesions have irregular shapes, some of them were regular, similarly most of the benign lesions have regular shapes but some of them were irregular. Although we could collect equal number of examples from both classes randomly, we have manually selected images that have different shape characteristics to compose more robust model for training and testing. Our dataset is made available publicly [23].

The images are 224 × 224 pixels and are labeled as benign and melanoma with binary specifiers. 2 examples are shown in Figures 2 (a) and (b).

Figure 2. (a) Benign skin image, (b) Melanoma lesion image

In this research, we focused on analyzing lesion images using Fourier Descriptors. Descriptors are the mathematical representation of the contour of an image using the Fourier transform. This technique captures the shape's essential features and can be used for shape analysis, recognition, and classification. Contour extracting operation is mostly tolerable to different acquisition conditions so from the frequency analysis point of view, model’s generalization ability is superior however color features are dependent on lighting. This is a big problem especially for deep learning models like CNN which processes the whole images and fully depends on their characteristics. On the other hand, ISIC is also working on establishing digital imaging standards. ISIC promotes the adoption of the Digital Imaging and Communications in Medicine (DICOM) standard for dermatology. DICOM, widely used in medical imaging fields such as radiology, cardiology, and oncology, enables standardized metadata, image formats, and device interoperability.

Figure 3. An unclear lesion image

Since contour plays the most important role, images must be cleaned to have correct form of contour. However, the lesion images obtained from ISIC database have hair, ruler shape and dark regions at the corners as shown in Figure 3.

Ruler and dark regions at the corners mostly were outside of the lesion and easy to clean but hairs generally cover lesion regions causing different contour formation therefore they must be cleaned before they are processed. To clean hairs, black hat filter method is applied [24]. This filter is used to remove dark regions over light regions. The background of images has always light color shines while hair has black therefore filter worked without a problem. A Python library function [25] has been used to apply black hat filter and has been executed for all the images which contain hairs on lesion.

3.2 Feature extraction

Feature extraction plays an important role in the overall success of the detection of cancer. In this section, first we propose a novel feature from FD’s energy and later add more selected discriminative features extracted from the shape of the lesions.

3.2.1 Fourier descriptor analysis

Fourier descriptors are derived by applying the Fourier transform to an image boundary, typically represented by its contour [26]. The resulting Fourier coefficients, known as Fourier descriptors (FD), defines the contour shape in the frequency domain. Low-frequency components describe the overall structure, while higher frequencies capture finer details of the shape. FDs are being calculated in 3 steps as:

Step 1: Contour extraction

Step 2: Expressing contour points by complex numbers

Step 3: Taking the Fourier transform from these complex numbers.

To extract contour images automatically, Sezgin and Taşaltı́n [27, 28] have implemented several effective methods defining thresholds and providing quantitative performance evaluation. Among those, Otsu’s method proposed by Otsu [29] has been used. It is a global thresholding technique that automatically selects the threshold value by minimizing the intra-class variance (or equivalently maximizing the between-class variance) between foreground and background. Firstly, image has been converted to grayscale with intensity levels [0, L-1] where L is 256. The goal here is find a threshold t that separates the image into two classes as C0 and C1:

C0: pixels with intensities [0, t] (background).

C1: pixels with intensities [t +1, L−1](foreground).

Then we compute the normalized histogram of the image using:

$p(i)=\frac{n_i}{N}$              (1)

where, n_iis the number of pixels with intensity i, and N is the total pixels.

For a candidate threshold t, class probabilities:

$\omega_0(t)=\sum_{i=0}^t p(i)$                  (2)

and

$\omega_1(t)=\sum_{i=t+1}^{L-1} p(i)$                (3)

Class means:

$\mu_0(t)=\frac{1}{\omega_0(t)} \sum_{i=0}^t i p(i)$                   (4)

$\mu_1(t)=\frac{1}{\omega_1(t)} \sum_{i=t+1}^{L-1} i p(i)$                   (5)

Otsu [29] defines the between-class variance as:

$\sigma_b^2(t)=\omega_0(t) \omega_1(t)\left(\mu_0(t)-\mu_1(t)\right)^2$                 (6)

The value of t that maximizes $\sigma_b^2(t)$ gives the optimal threhold value t* for the contour:

$t^*=\arg \max _t \sigma_b^2(t)$              (7)

Otsu [29] tries every possible threshold. The threshold that makes the groups most distinct (highest between-class variance) is chosen.

To express the contour points by complex numbers, let the contour of the shape be represented by a set of N points P₀, P₁, P₂, …… P_N where each point P_n = (x_n, y_n) is a pair of coordinates. We treat these coordinates as complex numbers:

$z_n=x_n+i y_n$                  (8)

where, i is imaginary unit. By applying the discrete Fourier transform (DFT) to this sequence of complex numbers, we obtain the Fourier descriptors Z_k as:

$Z_k=\frac{1}{N} \sum_{n=0}^{N-1} Z_n e^{-2 \pi i k n / N}$               (9)

where, k = 0, 1, ....., N-1.

The Fourier descriptors Z_k represent the frequency components of the shape, with lower-frequency descriptors capturing the overall shape and higher-frequency descriptors capturing finer details or irregularities. To make the Fourier descriptors invariant to translation, scaling, and rotation, we normalize them as follows:

We set Z₀ = 0 (the DC component) to remove the effect of translation. Then we normalize the descriptors by dividing them by |Z₁|, the magnitude of the first non-zero descriptor. We used the magnitudes |Z_k| of the descriptors, since the phase is affected by rotation. In this study, we have calculated 100 FDs using the explained method for 1000 images.

Irregularities in a shape’s contour are often captured by the higher-order Fourier descriptors. By analyzing the frequency domain energy distribution of the Fourier descriptors, we can identify whether a shape is more regular or irregular. The energy E of the Fourier descriptors is computed as:

$E=\sum_{k=1}^{N-1}\left|Z_k\right|^2$                     (10)

When we plot the energy of the descriptors in the frequency domain, we observe that the energy will be concentrated in the lower frequencies that do not differ much between benign and melanoma lesions; however for melanoma lesions which have irregular shapes, the energy has still strong components in the high frequencies as shown in Figures 4 (a) and (b).

Therefore, energy at higher frequencies can be used as discriminative property. Energy distribution shows exponential decay in relation to lower frequencies for both benign and melanoma lesions while more like linear at higher frequencies.

(a)

(b)

Figure 4. Energy spectrum of (a) benign (b) melanoma

3.2.2 Frequency change-point evaluation

To find the frequency change-point (threshold) from exponential decay to linear, we employed a dimensionality-reduced probabilistic modeling approach to detect change-points in multivariate time series data. First, Principal Component Analysis (PCA) [30] was applied to the dataset (100 FD series for 1000 image data set) to extract the first principal component, representing the dominant shared temporal trend. 1st principal component is taken for PCA. The related parameter is n_component and selected as 1. The principal components which explains 80% of the variance in the dataset were extracted and treated as a univariate representation of the collective behavior of all signals. This step reduces noise and aligns the temporal dynamics of the series into a single trajectory suitable for sequential modeling.

This univariate signal was then modeled using a Gaussian Hidden Markov Model (HMM) [31] with two latent states. The HMM captures changes in the statistical properties of the signal over time by assuming a latent Markov process and Gaussian emission distributions. The model was trained using the Expectation-Maximization (EM) [32] algorithm to maximize the likelihood of the observed sequence. The most probable sequence of hidden states was inferred using the Viterbi algorithm [33]. The change-point was defined as the first point at which the inferred hidden state transitioned. This method enables robust detection of regime changes in high-dimensional temporal data by leveraging both variance-based dimensionality reduction and probabilistic state modeling.

A Python program has been developed using hmmlearn library [34] to find the change-point. There were 4 parameters for HMM as model, n_components, covariance_type and algorithm. The model was selected as GaussianHMM, which uses HMM with Gaussian emissions. Since 2 latent states were selected, the parameter n_components was set to 2 which represents number of states. Covariance_type was set to “diag” which means each state uses a diagonal covariance matrix. Viterbi algorithm finds the most likely sequence of states, so algorithm parameter was set to “viterbi”. The change-point was identified as the first temporal index at which the predicted hidden state changed, indicating a statistically significant transition in the underlying generative process of the data. This transition corresponds to a qualitative change in the regime of the signal a shift from exponential decay to linear behavior. The change-point occurred at around 40 Hz, as shown in Figure 5.

Figure 5. Change-point detected by PCA + HMM

3.2.3 Novel feature extraction

After the change point, energy decay shows linear characteristic. Given that energy |Z_k|² changes as the frequency index k increases, we propose modeling the squared magnitudes of FDs as a function of the frequency index:

$\left|Z_k\right|^2=\beta_0+\beta_1 k+\epsilon_k \quad k=1,2, \ldots \ldots \ldots . N-1$               (11)

Here, β₀ is the intercept, β₁ represents the slope (rate of energy decay), and ϵ_k is an error term. For regular shapes, the slope β₁ will be more negative, indicating rapid energy decay. For irregular shapes, β₁ will be closer to zero, reflecting a slower decay or more evenly distributed energy.

We now introduce an OLS [35] model to analyze the energy decay across Fourier descriptors after 40 Hz. The OLS approach given in Eq. (11) fits the regression model. In the model, we minimize the sum of squared residuals:

$\tilde{\beta}_0, \tilde{\beta}_1=\arg \min \sum_{k=1}^{N-1}\left(\left|Z_k\right|^2-\left(\beta_0+\beta_1 k\right)\right)^2$                   (12)

The parameters $\tilde{\beta}_0$ and $\tilde{\beta}_1$ are estimated using standard OLS techniques by a Python program. Model parameter was selected as Linear Regression from Sklearn library. The key aspect of this model is the slope, $\tilde{\beta}_1$, which provides a quantitative measure of energy concentration.

Albay and Kamaşak [14] have used Fourier Descriptors as well. They have used limited number of Fourier Descriptors mainly in low frequencies achieving 83% accuracy but as we have shown with the spectral analysis in Figure 4(b), irregular shapes (mainly melanoma lesions) have significant components in the higher frequencies. Instead of including all frequency components, by this method we have managed to decrease hundreds of Fourier descriptors to one feature named as slope.

3.2.4 Statistical performance of this novel feature

Following hypothesis is used for the testing of the slope $\tilde{\beta}_1$:

H0: $\tilde{\beta}_1$ = 0 (no energy decay, high-frequency behavior)

H1: $\tilde{\beta}_1$ < 0 (significant energy decay, low-frequency behavior)

The test statistic is:

$t=\frac{\tilde{\beta}_1}{S E\left(\tilde{\beta}_1\right)}$                    (13)

where, SE($\tilde{\beta}_1$) is the standard error of $\tilde{\beta}_1$. If the test statistic is sufficiently negative, we reject H0, indicating significant low-frequency behavior.

A 95% confidence interval for $\tilde{\beta}_1$is given by:

$\tilde{\beta}_1 \mp t_{\frac{\alpha}{2}, N-2} . S E\left(\tilde{\beta}_1\right)$               (14)

If the confidence interval contains only negative values, it confirms that energy is concentrated in low frequencies, indicating a regular shape.

The goodness of fit is measured by the R² statistic, which captures the proportion of variance in |Z_k|² explained by the frequency index k:

$R^2=1-\frac{\sum_{k=1}^{N-1}\left(\left|Z_k\right|^2-\left(\widetilde{\beta}_0+\widetilde{\beta}_1 k\right)\right)^2}{\sum_{k=1}^{N-1}\left(\left|Z_k\right|^2-\left|\widetilde{Z}_k\right|^2\right)^2}$                 (15)

High R² indicates that the OLS model provides a good fit for the energy decay.

The fitted slope $\tilde{\beta}_1$ is the main discriminator for distinguishing between low and high-frequency behavior:

$\tilde{\beta}_1$ < 0 (dominant low-frequency components, regular shape)

$\tilde{\beta}_1$ ≈ 0 (balanced low and high frequencies, irregular shape)

Provided with the statistical testing, estimated $\tilde{\beta}_1$ is calculated for every image and selected as a primary feature for classification analysis.

3.3 The ANN model

An ANN is a computational framework inspired by the structure and functioning of biological network composed of many neurons. The flow of information through the network influences its structure, as the ANN adapts—or “learns”—based on the given inputs and resulting outputs. ANN consists of many simple, interconnected processors (neurons). It processes input signals by dynamically adjusting adaptive weights through a learning method known as backpropagation, which uses predefined training data and desired outputs. The network updates these weights by minimizing errors between predicted and actual outputs using an error function.

Our dataset was composed of 1000 224 × 224 pixel images from equally distributed between benign and melanoma cancer lesion. In addition to the extracted feature from FDs, average RGB (Red Green Blue) color values and shape parameters (major/minor diameters, eccentricity, compactness and perimeter) are calculated from images and provided as input to the ANN. We have used ANN with supervised manner by adding labels marking melanoma and benign skin lesions. Our moderate size model has 1 input layer, 2 hidden layers and 1 output layer created using the Tensorflow library in Python [36]. Early stopping is also used to improve the performance of the model. Model parameters are given in Table 1.

Model was performed with core i7 CPU @ 1.80 GHz processor with 8 Mb memory. Model execution has been completed in a few minutes.

Table 1. ANN layers and parameters

The performance of machine learning models is measured using the confusion matrix. As shown in Table 2, the error matrix consists of the following values: TP, which shows the number of true positives, FN, which shows the number of false negatives, FP, which shows the number of false positives, and TN, which shows the number of true negatives.

Table 2. The confusion matrix

Performance measures calculated from the confusion matrix are:

Accuracy $=\frac{(T P+T N)}{(T P+T N+F P+F N)}$                  (16)

Precision $=\frac{\mathrm{TP}}{(T P+F P)}$                 (17)

Recal $=\frac{\mathrm{TP}}{(T P+F N)}$               (18)

Specificity $=\frac{\mathrm{TN}}{(T N+F P)}$                 (19)

Simple and interpretable yet very effective pipeline for the detection of melanoma with a novel feature based on the spectral energy distribution of Fourier shape descriptors, called slope, supplemented by color and shape features of the lesion is proposed. The proposed efficiency has been approved by AI explainability tool SHAP. It has so balancing and separable power of melanoma and health lesion, it needs only a few features to be supported better accuracy. As far as we know, it has not been announced yet.

Slope feature efficiency leads to a small size ANN which doesn’t require more powerful computer resources for developing, even an applicable model to electronic devices with limited resources such as cell phones. The model accuracy is nearly as high as the best models obtained in literature in addition to such advantages it has. The reason is the feature we discovered, the spectral energy distribution of Fourier shape descriptors, what the scientist has been researching to find a shape feature on the lesion boundary they feel heuristically.

The contribution of current work has three important outputs to literature: a novel feature that brings powerful information to separate binary classes efficiently, the spectral energy distribution of Fourier shape descriptors; lightweight a machine learning model with minimized parameters; and high balanced class accuracy which is important in such health problems.

The only criticism for the features is the color properties which seem dominant in the classification with vision of SHAP analysis. However, color may have differences depending on the environmental conditions and devices used. Considering this instability, it would be better to further experiments with other data sets obtained from real hospital environments.

Even though lesion area edge features are very important, the dimension is also another issue. Three-dimensional vision brings the height of the lesion area, which is another important property that should be considered together with slope in further works.

Data is available online:

https://www.kaggle.com/datasets/cengizriva/benign

https://www.kaggle.com/datasets/cengizriva/melanoma

https://www.kaggle.com/datasets/cengizriva/feature-set

The codes which show the performance results of this article will be made available by the authors on request.