On the geometric modulation of skin lesion growth : a mathematical model for melanoma

Introduction: Early detection of suspicious skin lesions is critical to prevent skin malignancies, particularly the melanoma, which is the most dangerous form of human skin cancer. In the last decade, image processing techniques have been an increasingly important tool for early detection and mathematical models play a relevant role in mapping the progression of lesions. Methods: This work presents an algorithm to describe the evolution of the border of the skin lesion based on two main measurable markers: the symmetry and the geometric growth path of the lesion. The proposed methodology involves two dermoscopic images of the same melanocytic lesion obtained at different moments in time. By applying a mathematical model based on planar linear transformations, measurable parameters related to symmetry and growth are extracted. Results: With this information one may compare the actual evolution in the lesion with the outcomes from the geometric model. First, this method was tested on predefined images whose growth was controlled and the symmetry known which were used for validation. Then the methodology was tested in real dermoscopic melanoma images in which the parameters of the mathematical model revealed symmetry and growth rates consistent with a typical melanoma behavior. Conclusions: The method developed proved to show very accurate information about the target growth markers (variation on the growth along the border, the deformation and the symmetry of the lesion trough the time). All the results, validated by the expected phantom outputs, were similar to the ones on the real images.


Introduction
Cancer and related diseases remain a major issue to the scientific community.Due to the nature of the underlying process of carcinogenesis, several multicenter teams have been making efforts to further understand the process (at the molecular, cellular and tissue levels).The research in cancer provides evidence that an interdisciplinary approach brings empowerment to research, and this is an important contribution to what we know today about cancer and carcinogenesis.
Despite being the rarest among all skin cancers, the melanoma is one of the most deadly forms of the disease (Longo et al., 2012).Melanoma is a skin cancer derived from the melanocytes, cells that produce melanin, the skin coloration pigment.
The United States epidemiology of cancer laid skin cancer as the most common form of malignancy over the past three decades (Rogers et al., 2010;Stern, 2010).In particular the melanoma is also the most common form of cancer for young adults (25-29 years old) and the second most common form of cancer for young people (15-29 years old) (Bleyer et al., 2006).
Due to the high malignancy potential of melanoma, early detection of suspicious skin lesions is critical to prevent malignancy and to increase the treatment efficacy (Whiteman et al., 2011).
Dermoscopy is a non-invasive imaging technique used to obtain digital images on the skin surface by using a device known as dermoscope.It comprises a magnifying lens and a light source attached to a digital camera.This device allows the visualization of pigmented structures or vessels in the epidermis and superficial dermis (Zalaudek et al., 2010).Most dermoscopic structures are correlated with specific histopathological markers, therefore dermoscopy can be regarded as a link between clinical (macroscopic) and histopathological (microscopic) morphology (Zalaudek et al., 2006b).
Among all features that we can extract from images, as a consequence of the carcinogenesis evolution, the lesion growth is of particular relevance.Although this is not the only symptomatic feature, the growth and growth rate pattern of the lesion is one of the most valuable signs of malignancy, and can be studied using mathematical models that quantify and classifies the growth pattern.The combination of image processing techniques with mathematical models is a strategy with an important role to map the progression of lesions.
This work presents an algorithm based on linear transforms that evaluates the evolution of the border of the skin lesions based on two main measurable markers: the symmetry and the geometric growth path of the lesion.These markers are obtained through the relation between two segmented dermoscopic images of the same melanocytic lesion, obtained at different moments in time.
After this introduction of the problem, the remaining sections of this article are organised as follows: In the methodology and image setting section, methodological questions related to the image pre-processing are presented and the steps of the mathematical model are described.Then, in the section of results, the model is applied to a set of images (phantoms and real images).The work ends with the discussion and conclusions section, presenting the analysis of the results and concluding with the future perspectives in the context of the mathematical modelling of growth lesions in skin cancer, suggesting in perspective an individual signature model, in which growth and symmetry plays a pivotal role.

Methods
Dermoscopic images are obtained in the context of a clinical examination, and can be affected by some lack of time imposed by economic and human constrains.These and other aspects can affect the quality of the acquired image.Also the presence of hairs and air bubbles (when a fluid is applied in the interface between the dermoscope and the skin) requires a pre-processing step to assemble the image to the segmentation.The segmentation, being the first step in the image analysis, requires a robust methodology, since the remaining process depends on its accurate output.
The method used to segment images in this work was firstly proposed for dermoscopy in (Pereira et al., 2015), and uses the Multiscale Local Normalization (MLN) methodology.The use of a Gaussian kernel to the image in MNL overcomes the image perception dependence on the surroundings.
Matlab® software was used to implement the model presented in this work.It was used also in the construction of a set of binary images.
Real dermoscopic images used to test the methodology were obtained in the context of clinical practice by a dermatologist in a public hospital, with a previous informed consent to use these images for research purposes.

Image dataset
To test the robustness of the algorithm, a set of pairs of binary images, which were labeled as "phantoms" were created.This set of phantoms represents different shapes for image segmentation and different patterns of possible growth.Different scenarios of growth and growth rate pattern (regularity) and also symmetry are simulated in these images, covering a wide range of possibilities associated with these two morphological markers used in melanocytic image processing.The set of phantoms used is shown in Figure 1 in which it is possible to observe three different patterns of growth simulation: (a) growth preserving the shape; (b) growth maintaining shape but not proportional in every direction with respect to a given centroid; (c) and (d) growth changing shape and not proportional in every direction relative to a given centroid.

Image registration
According to what has been stated above, two images of the same lesion, obtained in two separate moments, were used to establish a marker of growth.It cannot be guaranteed that both images are captured exactly in the same conditions, i.e., the angle of incidence, light scattering and skin properties, among others.These different conditions may induce errors in the quantification of growth, caused by non-aligned pairs of images.In order to minimize these errors, a process of image registration was performed.
In this step the goal is to obtain the best alignment between both images.Hence, both images are aligned by their centroids.Among all possible senses for best alignment (mathematical, morphological, biological), taking into account the kind of images and the context of its acquisition, the global shape marker is used to achieve this matching, as described in the following.
The adjustment is performed by rotating the image taken in time one, I t 1 , over the image taken in time zero, I t 0 , and simultaneously calculating the standard deviation between them, defined by ( ) ( ) where ( ) x y and ( ) x y are respectively the points of the border of I t 0 and I t 1 obtained by the segmentation process.
The best adjustment is obtained by minimizing the standard deviation.In Figure 2 it is possible to observe an example of this adjustment in which (Equation 1) is minimized.

Border points
In order to optimize the process, an equally-spaced undersampling approach is performed on the border of the lesion I t 0 .This allows reducing the amount of data that the algorithm works, while preserving all the essential information, which basically is the lesion shape.For each of the selected points in I t 0 , there must be a correspondent one in I t 1 . Those matching points are selected to be the nearest points to the respective one in the border of I t 0 (see Figures 3a and b).In some cases, where the distance between border I t 1 candidates is abnormal (see Figure 3c), another approach is required.Normal/abnormal distances were classified according to previous simulation studies in other melanocytic images (Ferlay et al., 2013) from different image databases.In abnormal distances cases, the points will be selected as the intersection between the I t 0 border and a line that passes through the centroid ( ) x y and the point ( ) x y in I t 0 .

Global growth rate
The Global Growth Rate (GGR) defined in (2) reports information about the overall lesion growth behavior computing relative measures between the undersampled points.This rate is defined by where ( ) ( ) and 1, , for a sampling with M points where the points ( ) x y and ( ) x y are respectively points in the border of I t 0 , the corresponding points in the border of I t 1 and the centroid.
The different behavior of GGR is presented on Table 1A, where, T 1 = 0 in lesions with absence of growth, for all i = 1,…,M, is obtained directly from the GGR i expression, because if there is no growth then the numerator in (Equation 2) is zero for every i.To detect growth, the values of GGR i , i = 1,…,M, must be different than zero for one sample of these values.
The deviation of , denoted by ∆GGR, was also computed in order to verify if the lesion grows globally the same distance along each point on the border, and is given by (Equation 3).
In Table 1B, T 2 = 0.2 is obtained by experimental results between predicted and labeled phantom behavior.For every ∆GGR below this number is assumed to represent a lesion with a non-significant variation of the GGR values, therefore being labeled as having the same growth along the border.Otherwise, they are labeled as lesions with variation on the growth along a sample of points of the border.

Local growth rate
The Local Growth Rate (LGR) presented in (Equation 4) gives information on how the lesion grows locally and is given by ( ) ( ) If the lesion has a regular local growth behavior then it indicates absence of shape alterations.In other words, growth has been also proportional through time in every direction in relation to local dimensions.Otherwise, on the presence of irregular local growth, the lesion suffered deformation.If the deformation is small, it is not possible to detect alterations with naked eye, as presented in the example of Figures 5A and B.
This deformation is quantified by the deviation of , denoted by ∆LGR, which is defined in (Equation 5).
In Table 1C, T 3 = 0.2 is obtained by experimental results between predicted and labeled phantom behavior.For every ∆LGR bellow this number is assumed to have a non-significant variation of the LGR values, therefore being labeled as a lesion with maintained proportion.Otherwise, they are labeled as lesions with deformation., and (c) using radial lines as reference.Those matching points are selected to be see Figure 3a and b.In some cases, where the distance between border I t 1 candidates is abnormal (see Figure 3c), another approach is required.

Circularity
Circularity values defined in (Equation 6) inform about the circular type of the lesion in time zero, I t 0 and it is denoted by S and defined by the statistical coefficient of variation (CV) of elements in the set of coefficients between local and global growth, CLG.

S CV CLG CLG
where CLG denotes the average of elements of CLG, s CLG denotes the standard deviation of the elementsCLG, with CLG defined by , for LGR CLG GGR = .
As mention, circularity data informs about the circular type of the shape of the lesion in I t 0 and cannot have negative values.If GGR and LGR curves results have similar behavior then the lesion in I t 0 must have a circular type shape as both local and global distances produce the same results.In fact, the elements in the set CLG will mostly be 1.Therefore the value of S is closer to 0, which means there is no variability in the data in CLG.It is also possible to check information about the regularity of the shape of the lesion I t 1 , if in expressions (Equation 2) and (Equation 4) we replace d i0 by ( ) ( ) . In this case, we analyze the circularity shape of the lesion I t 1 .In Table 1 (D), T 4 = 0.1 is obtained by experimental results for predicted and labeled phantom behavior.Every lesion with circularity below this number is assumed to have a significant relation between LGR and GGR, therefore it is labeled as a circular type shape lesion.Otherwise, it is labeled as a non-circular shape lesion.

Symmetry coefficients
Symmetry is widely used to access the potential evolution from naevus to melanoma in dermoscopy skin evaluations.In order to access the lesion symmetry as an overall marker, the Angular Coefficient of Symmetry (ACS) is computed for each pair of images (I t 0 and I t 1 ).To calculate this coefficient an orthonormal reference system is considered whose origin coincides with the centroid of the image.Next it is obtained a symmetry linear transformation over the y-axis, as defined in (Equation 7).
where (x, y) are lesion points in the 1 st and 4 th quadrants.
Next an overlapping between ( , ' x y ′ ) and the lesion points on the 2 nd and 3 rd quadrants is obtained.Since the image is binary, after these overlapping, the final result will display 0, 1 or 2 for pixel values.Pixels with 0 indicate non-lesion with non-lesion overlapping, pixels with 1 indicates lesion with non-lesion overlapping and for those pixels with value 2, this indicates lesion with lesion overlapping.The relevant information is given by the relation between the frequencies of the pixel value 2 over 1.The higher the relative presence of 2, the higher the symmetry response, which is defined in (Equation 8).where z j = total number of values j in the 2 nd and 3 rd quadrants, { } 1, 2 j∈ .The symmetry over the x-axis is obtained using a linear transformation defined in (Equation 9) by where (x, y) are lesion points in the 1 st and 2 nd quadrants.
In the following an overlapping between ( , ' x y ′ ) and the lesion points on the 3 rd and 4 th quadrants is obtained.As in the former case a coefficient of symmetry SC x is computed and defined in (Equation 10).where z j = total number of values j in the 3 rd and 4 th quadrants, { } ACS can be used for quantitative and qualitative proposes, depending on the goal.For a qualitative usage, a threshold value is set in order to define the lesion as either symmetric or asymmetric, using the decision Table 1E.Otherwise, if the intent is to obtain a more discriminant measure, the value of ACS for each lesion may be used.
In Table 1E, T 5 = 0.8 is obtained from experimental results between predicted and labeled phantom behavior.Every lesion with ACS above this number is assumed to have a significant symmetry along their rotation, therefore it is labeled as a symmetric lesion.Otherwise, it is labeled an asymmetric lesion.
The deviation of ACS, denoted by ∆ACS, and given by (Equation 12) was also computed in order to correlate this results with the ones obtained for the regularity of the shape in circularity.If the deviation of ACS is small, that means that the ACS along the lesion is fairly constant, therefore the lesion must have circular shape with possible small variations along the border.Otherwise, it means the lesion has a non-circular shape.So when facing lesions with these characteristics, the decision should match the obtained previously on the circularity.
In Table 1F, T 6 = 0.2 is obtained by experimental results between predicted and labeled phantom behavior.Every lesion with ∆ACS below this number is assumed to have a non-significant variation of their symmetry; therefore it is labeled as a circular type shape lesion.Otherwise it is labeled as non-circular type shape lesion.

Results
In this section we discuss the results obtained using the markers defined above.First, the performance is tested in phantoms to validate their response, through the variation of growth and shape parameters.Then the procedure is tested on real images.
In Figure 4 there are two examples of lesions with absence of growth.There is a circular type shape phantom lesion in Figure 4a and a non-circular type shape phantom lesion in Figure 4b.This is just a preliminary step in order to validate the response of the proposed method on the absence of growth.
There is no growth when the GGR equals zero in every direction from their centroid, that is, when the GGR mean is 0 and the corresponding ΔGGR is also 0. As expected in this case, GGR mean and ΔGGR values from the two examples in Figure 4 are both zero.
Figure 5A presents the same phantom lesions used in Figure 4, but in this case both lesions show the same growth along the border.It can be observed that GGR in both images is quite regular, and this is corroborated by ΔGGR which in both cases is less than 0.2.Although naked eye observation shows that both lesions maintained their proportions, a careful examination of ΔLGR value shows that this behavior is only found for Figure 5A-a.In fact, in Figure 5A-b the lesion has suffered deformation, this is confirmed by their ΔLGR higher than 0.2.The circularity values confirm what is visually suggested, that is the image I t 0 has a circular type shape in Figure 5A-a and a non-circular type shape in Figure 5A-b.
In Figure 5B there are two examples of lesions with abnormal cases of growth.In Figure 5B-a a non-circular type shape lesion grows into a circular one, and in Figure 5B-b a circular type shape lesion grows into a non-circular form.It can be observed that GGR and LGR behavior in both images is quite irregular, which is confirmed by their ΔGGR and ΔLGR values above 0.2.This means that their absolute growth is irregular and their shape suffered deformation through time.Also in this case, the circularity values confirm what is visually suggested, that is image I t 0 is a non-circular type shape in Figure 5B-a and a circular type shape in Figure 5B-b.
In Figure 5C there are two cases of real lesions that evolved over time, obtained from two patients in a routine dermatologist examination.It is observed that both images have variations on their growth along the border, as both ΔGGR is higher than 0.2.
Regarding the local growth, and analyzing ΔLGR, both these lesions have values above 0.2 showing a presence of deformation.Furthermore, it is clear that in Figure 5C-b we are in the presence of a substantial deformation.
Finally, in what concerns circularity of the lesion we show two cases in Figure 5C.In Figure 5C-a, even though we have the presence of small variations on the border, the shape is circular as confirmed by circularity <0.2.In Figure 5C-b we are in the presence of a non-circular type shape.These two results are consistent with the visual assessment.
For global symmetry, this coefficient is close to 1 when the image is symmetrical or almost symmetrical (referred to the orthonormal reference system considered).As regards the phantom evaluation, notice that the values of ACS never exceeded 1.So when the average values of ACS, represented by ACS t 0 and ACS t 1 , are above 0.8 (Table 1E), the image is assumed to be symmetrical, regardless of the angle of rotation of the orthonormal reference system considered.This can be observed in Figure 6A-a for I t 0 and I t 1 , in Figure 6B-a for I t 1 and in Figure 6B-b for I t 0 .In the remaining phantom images, the lesions lack symmetry and the corresponding average values of the ACS decrease.
It is also calculated the deviation of the ACS in each case, denoted by ∆ACS t 0 and ∆ACS t 1 .Notice that, if the symmetry of the lesion is independent of the angle of rotation of the orthonormal reference system considered, these coefficients are closer to zero, and therefore the lesion has a circular type shape.Following Table 1F, this can be observed in Figure 6A-a for I t 0 and I t 1 , in Figure 6B-a for I t 1 and in Figure 6B-b for I t 0 . If the deviation of the ACS is higher than 0.2, as verified in the rest of the phantoms, this means that there are some specific positions (angles) of the orthonormal reference system considered for which the lesion is more symmetric, and therefore it is not of a circular type.
In Figure 7 it is considered the case of two real lesions that evolved through time.When comparing the two images, it is possible to observe that Figure 7a is visually symmetrical while Figure 7b is not.That is reflected on their ACS t i , showing in Figure 7a a symmetric lesion and in Figure 7b lack of symmetry.
Finally, ∆ACS t i is also consistent with the decision for circularity, as in Figure 7a it is presented an image with a deviation smaller than 0.2, therefore is assumed to be of a circular type shape, and in Figure 7b we have the opposite.

Discussion
Early diagnosis is important for the prognosis of melanocytic skin lesions and the measurements obtained through the proposed mathematical model increases the discriminating capability of the classical diagnosis when applied at the follow-up of the lesions.The ability to quantify the growth pattern of a certain lesion is important when accessing growing and growth rate.The method developed proved to show very accurate information about the target growth markers (variation on the growth along the border, the deformation and the symmetry of the lesion trough the time).All the results, validated by the expected phantom outputs, were similar to the ones on the real images.Also the measurements of the circularity and the ∆ACS, which were derivate from the previous markers and both used to quantify the circular type of the shape of the lesions, not only have coincident qualitative results, but also have the expected response for all the images.
The assessment of these quantitative markers can be used to increase the efficacy of automatic classification using machine learning algorithms.
The use of growing factors and symmetry coefficients joint with texture analysis features is now being considered and worked in our research group.We have to acknowledge to Juliana Baptista MD, Dermatologist for her contribution in image acquisition.

Figure 1 .
Figure 1.Phantom examples for growth simulation (a) preserving the shape; (b) maintaining shape but not proportional in every direction; (c) and (d) changing shape.

Figure 2 .
Figure 2. Example of image registration.Left: Images I t 0 and I t 1 before the alignment.Right: same image after the alignment process.

Figure 3 .
Figure 3. Point detection approach in the border of image I t 1 .(a) and (b) by selecting the nearest points of I t 0

Figure 4 .
Figure 4. Phantom lesions with absence of growth.(a) circular phantom; and (b) non-circular phantom.

Figure 5 .
Figure 5. Growth Factors.A: Growth factors for phantom lesions with regular absolute growth -(a) regular shaped and (b) irregular shaped; B: Growth factors for phantom lesions with irregular absolute growth -(a) evolves from an irregular shape to one regular and (b) evolves from an regular shape to one irregular; C: Growth factor for real lesions.

Figure 6 .
Figure 6.Symmetry factors.A: Symmetry for phantom lesions with regular absolute growth -(a) regular shaped and (b) irregular; B: Symmetry for phantom lesions with irregular absolute growth -(a) evolves from an irregular shape to one regular and (b) evolves from a regular shape to one irregular.

Figure 7 .
Figure 7. Symmetry for real lesions.(a) example using a real image of circular type shape; and (b) example using a non-circular real image.

Table 1 .
Decision table for the different markers.

.
The symmetry coefficients were obtained through decile rotations of the image to complete 360 degrees.In each iteration of the process results for angles