Hough Transform Sensitivity Factor Calculation Model Applied to the Analysis of Acne Vulgaris Skin Lesions

Abstract

The number of infectious spots or pathological structures recorded on dermatological images is a tool to aid in the diagnosis and monitoring of disease progression. Dermatological images for the detection and monitoring of the evolution of acne infections are evaluated globally, comparing whether the increase or decrease in infectious lesions appearing on an image is significant. This evaluation method is only indicative since its accuracy is low. The accuracy problem could be improved by an exact count of the number of structures and spots appearing on the image. The mathematical function circular Hough transform (CHT) function implemented in MATLAB is here applied to develop a procedure for counting these structures. CHT has been used in the recognition of benign and distorted red blood cells, in the detection of pellet sizes in industrial processes and in the automated detection and morphological characterization of breast tumor masses from infrared images, as well as for the detection of brain aneurysms and use in magnetic resonance imaging. The sensitivity factor is one of the many parameters required to feed the CHT algorithm. Its choice is unclear as there is no proper methodology to select an optimum value suitable for each image. In this work, a procedure for determining the optimal value of the sensitivity factor is proposed The approach is validated by comparison with the results of the manual counting of the points (ground truth).

1. Introduction

To date, there is no standardized way of monitoring the evolution of the pathology caused by acne vulgaris infection. Although the best way of monitoring is to count the number of infected pores appearing on the skin, this is not usually performed in dermatological clinics. In addition, it has been shown that infected spots have a different fluorescence than pores that appear to be infected but are not.

The circular Hough transform function (CHT) is one of the most robust spot detection and counting algorithms available, but for its correct implementation it is necessary to choose a sensitivity factor, which until now has been chosen randomly or empirically.

In the present study, an automated way of choosing the optimal sensitivity factor for the algorithm and for each image to be processed is developed. All this is studied on fluorescence images of the disease.

As a result, we obtain an algorithm that provides us with the optimal sensitivity factor for counting the number of infected points that appear in fluorescence images processed with the CHT function.

Until now, the implementation of the Hough transform function in MATLAB had as a weak point the choice of the optimal sensitivity factor. This problem is solved in the present work and will allow the algorithm to be applied accurately by the entire scientific community.

When monitoring the course of a pathological processes in dermatology, such as acne infections [1], it is common to record images of the control [2].

The structures of interest that are analyzed in the dermatological images studied in this article have circular or quasi-circular shapes. It is therefore interesting to study the application of CHT in the processing of monochromatic dermatological images.

The CHT function is implemented in the MATLAB software, within the Image Processing Toolbox. This function has different input arguments such as radius or maximum and minimum diameter of the structures to be detected, background shape (brightness–blackness), etc.

Of all the input arguments, and as described below, one of the most difficult parameters to choose for the application of this function is the sensitivity factor S. This factor has values in the range (0–1). The closer we get to the value 1, the more faint and partially obscured circles the algorithm can detect, but at the same time the risk of false detection increases. With the same sensitivity factor value, we do not obtain the same results in two different images and, therefore, the choice of an optimal value for each image to be processed is important. Thus far, this value has been chosen empirically [3], through trial and error stack overflow.

As this way of choosing has proved to be inefficient and imprecise, and as the need to determine the optimum value for processing each image is evident, a methodology is proposed that allows a standardized and automated way of finding an optimum value for each image to be processed and analyzed. The possibility of quickly and automatically quantifying the number of structures that appear in medical images has been studied in works such as [4], where different segmentation algorithms are proposed. In [5], an algorithm for the recognition and quantification of cells present in different histological slices is proposed. In [6] it is used circle Hough transform to analyze immune responses by quantifying the amount of cytokine secreted by cells. ELISpot and FluoroSpot images were analyzed. [2] used the circle Hough transform for cell detection in blood samples. Ref. [7] used CHT for femoral head detection in perthes MR slices. Ref. [8] used CHT for the detection of abdominal aortic aneurysm automatically.

Although segmentation of medical images is addressed in many studies as [9,10,11,12,13,14,15], and segmentation of elements in inhomogeneous images is achieving very good results as in [16,17,18,19,20], quantification of the exact number of points in a dermatological image is a problem that has not been solved efficiently. In this paper, CHT is applied to count the number of points in dermatological images and the algorithm is implemented in MATLAB. The main problem with this algorithm is the choice of the optimal sensitivity factor for each image. From a given sensitivity factor, the detection of dots increases exponentially, the latter being the optimal factor for achieving a correct count.

As mentioned above, for each image we obtain several sensitivity factors, one of which gives us the highest precision in the point count. For each image, consecutive sensitivity factors are evaluated and for each sensitivity factor a number of points is obtained. The behavior of the data is as follows: from one sensitivity factor onwards, the number of points counted increases exponentially. This last sensitivity factor (before the exponential increase) is the optimum for the application of the algorithm.

The values of the sensitivity factor and the corresponding number of points obtained have been plotted on a cartesian graph (X,Y). A sharp change in the slope of the graph was observed when the number of points started to grow exponentially.

To choose the optimal sensitivity factor, a statistical methodology presented in studies such as [21] has been used. The behavior of the slope change can be modelled by choosing the mathematical equation according to the visual aspect of the XY plot or by knowing the expected behavior. In this study, the optimal point of the sensitivity factor will be sought by analyzing a change in slope. An algorithm for choosing the optimal sensitivity factor is proposed.

2. Background

2.1. Basic Idea of Circular Hough Transform Function (CHT)

Hough in 1962 proposed the representation of a line in space (r, Ө), so that straight lines could be detected in an image. Subsequently, this transform was extended to the detection of regular shapes. Practical cases of application of the Hough transform (HT) for linear elements can be found in line detection in low-level vision systems [22] and in the detection of straight lines in images with optimized algorithms [23].

The circular Hough transform (CHT) function is derived from the Hough transform (HT) for linear object detection. CHT is based on the same principle: a point (x,y) is represented by a curve in space (r, Ө) [24].

The equation of a circle answers to

(x − a)² + (y − b)² = r

(1)

Now, each centroid of circle (a, b) will be represented on a surface in the space (a, b, r,) where r is the radius of the circle.

The circle is thus defined by three parameters and, consequently, the parameter space will have three dimensions. A circle in image-gen-space will be a point in parameter space [25].

The parameter space is defined in each case by the shape of the object of interest to be detected in an image [26].

For each contour pixel there is a family of circumferences passing through that pixel and given by the parametric representation of circles

x = a + r* cos (t)

(2)

y = b + r* sin (t)

(3)

The angle t is evaluated in 0–360°.

The algorithm is based on a voting accumulation system [27]. Each contour pixel votes for all possible circumferences passing through that pixel. If a point appears in the parameter space that has accumulated many votes, the parameters of that point are interpreted as corresponding to the circle that passes through a large number of contour points [25].

Thanks to the cumulative voting system (maximum likelihood election) the CHT function is robust to noise and to the lack of data for some incomplete circles.

2.2. CHT in Matlab

CHT is implemented in the MATLAB software, ToolBox Image Segmenter in the function imfindcircles so that the parameters to be entered are the minimum and maximum diameter of the circular objects to be identified and a sensitivity factor, among others.

The sensitivity factor represents the minimum quality threshold of elements classified as circles. The value of the sensitivity factor must be in the range (0–1), the default being 0.85. With the same value of this factor applied to different images, the results obtained differ completely. By increasing the value of the factor, it will be possible to classify elements less defined in contour and elements less defined in contrast of intensity values with respect to their adjacent pixels [28]. As mentioned above, as the sensitivity factor value approaches 1, more faint and partially obscured circles can be detected with the algorithm, but the risk of false detection of elements also increases.

The main problem encountered is that there is no mathematical formula or method that allows a choice of the factor value that corresponds to an optimal performance of the algorithm. Thus far, the optimal value of the sensitivity factor has only been found empirically [3]. There is no methodology to automatically select an optimal value for each image to be processed. Therefore, we propose here a procedure for determining this value. The approach will be validated by comparison with the results of the manual counting of the points (groundtruth).

3. Material and Methods

MATLAB software version 2018b licensed by the Universitat Politècnica de València was used to process the images and implement the methodology for the automatic selection of the optimal value of the sensitivity factor. Hardware was PC Lenovo L380 Yoga.

The processed images were in raw format and the color model was monochromatic. The images were taken by illuminating the skin with ultra-violet light that activated the fluorescence of the lesions [29]. The fluorescence of infected lesions is caused by the bacterium Propionibacterium Acnes [30].

The images were taken with a SCHOTT RG610 high-pass optical filter placed in front of the sensor. The sensor was an IphoneX camera. As we only recorded the signal in the red range, the images were taken in monochromatic mode.

The groundtruth images used by the manual evaluators to count the points were printed on a Ricoh Aficio MP2550 laser printer with monochrome printing mode and 600 × 600 dots per inch resolution.

The manual dot counting was carried out by five different evaluators from the Centre for Research in Graphic Technologies of the Universitat Politècnica de València.

The methodology is described with reference to one of the images analyzed, Figure 1.

A small algorithm has been created that works in the form of a loop that for each image to be processed, applies consecutive sensitivity factors with an increment of 0.01, and for each one of these it registers the number of points detected in the binary image that the function automatically creates to carry out the count. The final output of the algorithm is a file with the number of points associated with each sensitivity factor for which the circular structure count has been performed. The algorithm will work for sensitivity factors between 0.81–0.82 up to a maximum of 0.96. Sensitivity values (S_factor) exceeding the limits of the range (0.81–0.96) have in no case (previous tests) given valid counts.

As mentioned, with the CHT algorithm implemented in MATLAB, the image will be processed by giving the sensitivity factor S values from 0.82 to 0.96 in 0.01 increments. As a result, a binary image (BW) will be obtained from which the number of identified points will be counted.

To count the number of elements identified by the algorithm, the function implemented in MATLAB bwlabeln is used, the input parameter of the function being the BW image obtained from imfindcircles.

The number of points detected associated with each sensitivity S_factor is recorded. With the different values of the number of points obtained, a new parameter (W) is calculated.

For each image, the algorithm is evaluated under different S_factors.

Table 1. Number of points detected for each S_factor.

Sfactor	Number of Points
0.82	403
0.83	403
0.84	471
0.85	605
0.86	605
0.87	657
0.88	708
0.89	787
0.9	877
0.91	945
0.92	966
0.93	1001
0.94	1087
0.95	1182
0.96	1279

shows the number of points for different S_factor. Then, for each S_factor that the algorithm evaluates, we obtain a number of points and a W. For each image evaluated, we obtain a set of data that are represented in the fitting model and from this model we choose the value of the S_factor under which the algorithm obtains a better count of the infected points.

For example, in Table 2 we can see how for each S_factor a number of points and a W are associated.

$$W=\frac{N{P}_i}{N{P}_{max}-N{P}_i}$$

(4)

NP_i = number of points for each sensitivity value.

NP_max = maximum number of points for maximum sensitivity value.

Table 2. Model fit for choice of S_factor in Trial 1.

Trial 1	SfactorR	0.81	0.82	0.83	0.84	0.85	0.86	0.87	0.88	0.89	0.9	0.91	0.92	0.93	0.94	0.95
	NP	-	403	471	549	605	657	708	754	787	877	945	966	1001	1087	1182
	W	-	0.460	0.583	0.752	0.898	1.056	1.240	1.436	1.600	2.182	2.829	3.086	3.601	5.661	12.190
	R2	-	-	94.10	95.26	96.46	97.62	98.64	99.39	99.56	98.16	96.06	94.74	-	-	-
	Signif.	-	-	X	X	X	X						X	-	-	-

X, any p-value ≥ 0.05 then the model is not valid; , all the p-values < 0.05 then the model is significant.

Table 2. Model fit for choice of S_factor in Trial 1.

Trial 1	SfactorR	0.81	0.82	0.83	0.84	0.85	0.86	0.87	0.88	0.89	0.9	0.91	0.92	0.93	0.94	0.95
	NP	-	403	471	549	605	657	708	754	787	877	945	966	1001	1087	1182
	W	-	0.460	0.583	0.752	0.898	1.056	1.240	1.436	1.600	2.182	2.829	3.086	3.601	5.661	12.190
	R2	-	-	94.10	95.26	96.46	97.62	98.64	99.39	99.56	98.16	96.06	94.74	-	-	-
	Signif.	-	-	X	X	X	X						X	-	-	-

X, any p-value ≥ 0.05 then the model is not valid; , all the p-values < 0.05 then the model is significant.

3.1. Statistical Analysis

Each W value will therefore be associated with a particular sensitivity factor S. The set of sensitivity factors S together with the respective associated W will be plotted as in Figure 3A As can be seen in Figure 2, a theoretical graphical representation with differentiable trend intervals is obtained. Thus, at the point on the graph where we find the beginning of an exponential interval, we will find the value of S that will be close to the optimum. Once the results of the count are obtained, the data will be fitted to a mathematical model with change of slope (SPLINE) so that the selection of the optimal point S will coincide with the X_R value that provides the best fit to the model in terms of R² and significance.

The theoretical slope change equation (SPLINE) is

Y = β₀ + β₁*X + β₂*D(X − X_R) + U

(5)

where:

Y = The dependent variable of the model.

β₀ = The independent coefficient. Represents the value of Y variable when the variable X is zero.

β₁ = The coefficient that accompanies the X variable. Represents the first slope of the model and measures the increase in the variable Y when the X increases one unit.

X = The independent variable of the model.

β₂ = The coefficient that accompanies the Dummy variable. Difference in the increase in the variable Y when one unit of the variable X is increased while X ≤ X_R with respect to when X > X_R.

D = Dummy variable that takes a value of 0 when X ≤ X_R, and a value of 1 when X > X_R.

X_R = The node from which the slope of the line will change.

U = The perturbation of the model. Represents the variables that do not appear in the model and that can affect the dependent variable Y.

Applied to the variables we are working with, this would be

W = β₀ + β₁*S_factor + β₂*D(S_factor − S_factorR) + U

(6)

where:

W = The dependent variable of the model.

S_factor = The independent variable of the model.

D = Dummy variable that takes a value of 0 when S_factor ≤ S_factorR, and a value of 1 when S_factor > S_factorR.

S_factorR = The node from which the slope of the line will change.

Table 2 shows the experimental data obtained for the trial 1 and the results for the fitting model. Every pair of R² and significance analysis were obtained using Equation (6), choosing each time the value of S_factorR that appears in the upper part of this pair. By means of the analysis of variance, the significance of the model and its parameters are checked, and by means of the coefficient of determination R², the fit of the model is evaluated. Subsequently, the results of S_factorR vs. R² are plotted and the value corresponding to the highest R² are assigned as the optimum point of S_factor, which in the graph appears at the highest point, as can be seen in Figure 3B.

3.2. Flowchart

In Figure 4 it is shown the global process for obtaining the optimal S_factor. For each image, the algorithm is evaluated under different S_factors. For each S_factor that the algorithm evaluates, we obtain a number of points and a W. For each image evaluated, we obtain a set of data that are represented in the fitting model and from this model we choose the value of the S_factor under which the algorithm obtains a better count of the infected points.

For example, in Table 2 we can see how for each S_factor a number of points and a W are associated.

4. Results and Discussion

Six different dermatological images shown in Figure 5 were analysed using the proposed methodology. Five different subjects were selected for manual counting of infectious spots on the images.

As we can see, we have produced a method that significantly improves and standardizes the process of choosing the S_factor. This improvement translates directly into:

-

An increase in the final precision in the counting of points by applying CHT.

-

Being an objective method, it allows the results of the algorithm to be compared with other point counting algorithms.

With the results obtained from the previous tests,

Table 3. Model fitting for choice of S_factor in the 6 patients.

Trial 1	SfactorR	-	0.82	0.83	0.84	0.85	0.86	0.87	0.88	0.89	0.9	0.91	0.92	0.93	0.94	0.95
	NP	-	403	471	549	605	657	708	754	787	877	945	966	1001	1087	1182
	W	-	0.460	0.583	0.752	0.898	1.056	1.240	1.436	1.600	2.182	2.829	3.086	3.601	5.661	12.190
	R2	-	-	94.10	95.26	96.46	97.62	98.64	99.39	99.56	98.16	96.06	94.74	-	-	-
	Signif.	-		X	X	X	X						X	-	-	-
Trial 2	SfactorR	0.81	0.82	0.83	0.84	0.85	0.86	0.87	0.88	0.89	0.9	0.91	-	-	-	-
	NP	187	205	220	236	256	285	320	378	471	662	985	-	-	-	-
	W	0.234	0.263	0.288	0.315	0.351	0.407	0.481	0.623	0.916	2.050	-	-	-	-	-
	R2	-	91.71	94.03	96.36	98.30	99.15	98.15	-	-	-	-	-	-	-	-
	Signif.	-	X	X	X				-	-	-	-	-	-	-	-
Trial 3	SfactorR	0.81	0.82	0.83	0.84	0.85	0.86	0.87	0.88	0.89	0.9	-	-	-	-	-
	NP	298	319	360	396	454	521	618	746	1015	1271	-	-	-	-	-
	W	0.306	0.335	0.395	0.453	0.556	0.695	0.946	1.421	3.965	-	-	-	-	-	-
	R2	-	87.09	90.64	94.23	97.12	98.81	96.71	-	-	-	-	-	-	-	-
	Signif.	-	X	X	X				-	-	-	-	-	-	-	-
Trial 4	SfactorR	-	0.82	0.83	0.84	0.85	0.86	0.87	0.88	0.89	0.9	0.91	0.92	0.93	-	-
	NP	-	47	65	108	140	167	198	253	323	430	571	688	859	-	-
	W	-	0.058	0.082	0.144	0.195	0.241	0.300	0.417	0.603	1.002	1.983	4.023	-	-	-
	R2	-	-	93.71	94.74	96.19	98.00	99.55	98.59	-	-	-	-	-	-	-
	Signif.	-	-	X	X	X				-	-	-	-	-	-	-
Trial 5	SfactorR	-	0.82	0.83	0.84	0.85	0.86	0.87	0.88	0.89	0.9	0.91	0.92	0.93	0.94	0.95
	NP	-	208	255	307	347	386	405	427	443	479	508	525	545	596	635
	W	-	0.423	0.573	0.781	0.983	1.229	1.373	1.564	1.724	2.167	2.646	3.000	3.516	5.731	9.769
	R2	-	-	96.42	96.91	97.42	97.98	98.72	99.42	99.91	99.54	98.60	97.63	-	-	-
	Signif.	-	-	X	X	X								-	-	-
Trial 6	SfactorR	-	0.82	0.83	0.84	0.85	0.86	0.87	0.88	0.89	0.9	0.91	0.92	0.93	0.94	0.95
	NP	-	276	300	320	340	355	374	399	423	448	473	497	541	582	618
	W	-	0.685	0.792	0.891	1.003	1.096	1.226	1.425	1.652	1.939	2.296	2.731	3.920	6.000	10.131
	R2	-	-	94.04	95.41	96.76	98.10	99.13	99.59	99.47	98.56	96.57	-	-	-	-
	Signif.	-	-	X	X	X							-	-	-	-

X, any p-value ≥ 0.05 then the model is not valid; , all the p-values < 0.05 then the model is significant.

and Figure 6 and Figure 7, together with the data obtained from a manual point count of the processed images (groundtruth), a correspondence factor has been obtained.

The number of points detected for each optimal S found has been compared with the number of points identified manually (obtained as an average of the manual count performed by five different evaluators under the same lighting conditions), the results of which are shown in

Table 4. Results for Sensitivity.

Trial	1	2	3	4	5	6
	804	263	379	166	442	405
	797	260	385	162	439	420
NP manual	802	268	389	169	440	417
	800	261	389	173	437	425
	808	259	392	163	447	431
$\overline{N{P}_m}$	802	262	387	167	441	420
Sfactor	0.89	0.85	0.84	0.86	0.89	0.89
NP for selected S	809	256	396	167	443	423
Sensitivity	0.992	0.976	0.976	0.998	0.995	0.992

.

The SEN factor has been obtained, which indicates the sensitivity of the methodology for the determination of S with the proposed methodology

$$SEN=1-Abs\left|1-\frac{NP for selected S}{NP manual}\right|$$

(7)

5. Conclusions

With the application of the proposed methodology, the problem presented by the empirical choice of the sensitivity factor for the application of CHT implemented in MATLAB has been solved.

Considering the sensitivity factors presented in the article, which have been selected based on the previous methodology, widely satisfactory results have been obtained in the counting of points (sensitivities of around 95%), as can be seen in the Table 4. Therefore, the proposed method is effective since it allows the determination of an optimum sensitivity factor automatically for each image to be processed. The method provides metric accuracy to the analysis of the images. Therefore, the exact number of points can be used for diagnosis instead of a global evaluation of them.

Further research is already ongoing to evaluate the robustness of the proposed methodology by means of a node similarity analysis between different graphs, as described by reference [31].

The possibility of finding the optimum value of the sensitivity factor for each processed image will allow the CHT function in MATLAB to be used for the study of images for which absolute point detection accuracy is required. An example of this is the work that has been initiated by the research centre in the field of detecting defects on ceramic tiles caused by excessive porosity. CHT is being applied because of the possibility of objectively obtaining an S_factor that ensures the exact count of the number of points.