Fast weighted least squares for detail and tone enhancement of medical images

Abstract

Background

The limitations in current medical imaging systems often result in poor visibility and insufficient detail. Such low-quality images hinder the decision-making of medical professionals. The weighted least squares (WLS) model provides a prevalent framework for detail and tone manipulation, but it is non-trivial to solve.

Method

We propose a fast solution to the WLS model. Instead of addressing a large linear system or approximating two-dimensional image processing with multiple one-dimensional procedures, we demonstrate that the model can be solved efficiently using an iterative algorithm. We leverage our method for the detail and tone enhancement of medical images.

Results

(1) For detail enhancement, on the three subsets of X-ray images, the spatial-spectral entropy-based quality (SSEQ) values are 21.15, 29.46, and 31.36; the convolutional neural networks image quality assessment (CIQA) values are 17.67, 25.60, and 26.01; on the computed tomography and magnetic resonance imaging images, the SSEQ and CIQA values are 14.01 and 16.93; it achieves peak signal-to-noise ratio (PSNR) ≥39.39 and structural similarity index (SSIM) ≥0.99; (2) for tone enhancement, the tone mapped image quality index, fidelity, and naturalness values are recorded at 0.7905, 0.9300, and 0.0090; it achieves PSNR ≥50.09 and SSIM ≥0.99. Notably, our method processes a 720P color image in 0.0725 seconds on a modern graphics processing unit.

Conclusion

The proposed method enhances the detail and tone of medical images. It addresses computational challenges inherent in traditional approaches. Its adaptability to diverse applications suggests a promising avenue for improving diagnostic accuracy and patient outcomes in clinical practice.

Keywords

Medical images detail and tone enhancement weighted least squares iterative algorithm

Introduction

Medical images such as X-rays, computer tomography (CT) scans, magnetic resonance images (MR), laser scanning confocal microscope (LSCM), and ultrasound are vital sources of detailed biological and anatomical tissue information crucial for clinical diagnosis and treatment by healthcare professionals.¹ Capturing conditions and technical limitations often lead to deficiencies in medical images, such as low contrast and poor visibility. Studies of Philip et al.² and Macgillivray et al.³ have indicated that a significant amount of medical images are not of sufficient quality for screening and diagnosis. Furthermore, the deficiencies in medical images would impede the subsequent image analysis tasks.

Driven by the aforementioned factors, there is a pressing need to overcome the deficiencies in medical images. The most popular traditional medical image enhancement methods are based on histogram equalization⁴ and its variants.^5–7 The histogram-based methods are efficient, but they fail to consider the spatial information and thus may be subjective to the discontinuity between neighboring regions. Another type of traditional medical image enhancement methods rely on the edge-preserving filter, such as the bilateral filter,⁸ guided image filter,⁹ and bilateral weighted least squares (WLS) filter.¹⁰ However, these filters may exhibit the halo or gradient reversal artifacts. The last type of medical image enhancement methods is based on neural networks, which have been widely applied to MRI and CT images.^11–13 The deep learning-based methods are effective, but they typically require a large amount of labeled data to train the model, which is time- and resource-intensive.

As a global image processing technique, the WLS model¹⁴ provides a prevalent framework for detail and tone manipulation, as it produces promising image quality and does not suffer from the halo and gradient reversal artifacts. However, due to the large sparse inhomogeneous Laplacian matrix in the linear system, the WLS model can be rather time-consuming to solve. To accelerate the image processing based on the WLS, researchers attempt to approximate the two-dimensional (2D) process as multiple one-dimensional (1D) processes.^15,16 The strategy certainly improves the efficiency, but it would also introduce the streak artifacts and sabotage the quality. The recent development of the deep learning technique provides the researchers with a new tool to explore the possibility of approximating the results of the WLS with neural networks.^17–20 However, these methods are either inflexible in tuning parameters or poor in processing quality.

The objective of this article is to develop an efficient and accurate solution to the WLS model for the detail and tone enhancement of medical images. Different from the histogram-based methods,^4–7 the proposed method considers the spatial information, thus not being prone to discontinuity between neighboring regions. Unlike the state-of-the-art image enhancement methods based on neural networks,^21–23 our method is free of the tedious labor used to define the vague labels of the images. Different from existing filter-based methods,^8–10 the proposed method is based on the WLS of the gradients, thus alleviating the halo and gradient reversal artifacts. The naive solution to the WLS model can be time-consuming, thus we propose a fast WLS solution. Unlike existing 1D approximate solutions,^15,16our method focuses on the solution in the 2D space, thus being immune to the streak artifacts. Compared to the learning-based approximate solutions,^17–20 our method supports parameter tuning, thus being more flexible.

The main contributions of this article are:

We decouple the gradient regularization from the weighted scheme. It breaks the foundation of the inhomogeneous Laplacian matrix that leads to high computation costs.

We propose a novel iterative solution to the WLS model, which is highly efficient.

Our method outperforms existing fast WLS methods in various quality evaluation metrics on the detail and tone enhancement of medical images.

The rest of this article is organized as follows. The “Literature review” section covers the works closely related to our method. The “Proposed method” section provides detailed descriptions of the proposed method. In the “Applications” section, we demonstrate how the proposed method can be applied to the detail and tone enhancement of the medical images. The experimental settings and results are presented in the “Experiment settings” and “Results” sections, respectively. The “Discussion” and “Limitations” sections are presented afterwards. Finally, we conclude the article in the “Conclusion and future work” section and present possible directions for future research.

Literature review

The WLS model¹⁴ is among the most prevalent frameworks for the tone and detail manipulation of images. It aims to minimize the sum of squared errors regularized by the weighted sum of the squared gradients, which is mathematically equivalent to the anisotropic diffusion model²⁴ with a nontrivial steady state. A key advantage of the model lies in that it substantially mitigates halo effects and gradient reversal artifacts. The suppression of these artifacts significantly improves the overall visual quality of enhanced images. However, despite its effectiveness, the practical implementation of the WLS model presents computational challenges. Specifically, the model requires the solution of a large linear system involving a sparse Laplacian matrix, which poses significant demands on both processing time and memory resources. The computational burden limits the scalability and applicability of the model in real-time or resource-constrained environments.

Numerous efforts have been dedicated to accelerating the solution of the WLS model, driven by the need to address its computational inefficiencies. Min et al.¹⁶ approximated the linear system defined over 2D space with multiple linear systems that are defined over 1D space, so that the large sparse coefficient matrix would become a more manageable tridiagonal matrix, which can be efficiently solved using established numerical techniques. Similarly, Liu et al.¹⁵ seek to approximate the 2D image processing task with multiple 1D ones. However, rather than constructing a four-connected or an eight-connected neighborhood, they adopt a zig-zag neighborhood construction scheme. This innovative strategy allows them to incorporate essential 2D neighborhood information while still capitalizing on the efficiency of 1D processing. Both of the two methods^16,15 achieve significant speedups, but the separated 1D processing would cause the streak artifacts.

In recent years, inspired by the booming advancement of deep learning techniques, researchers have started to investigate learning-based methods to approximate the results of the WLS model. A notable contribution in this domain is by Xu et al.,¹⁸ who propose a neural network architecture that learns the gradients of the output in a supervised manner. These learned gradients are subsequently utilized to reconstruct the images. Further contributions are illustrated by Chen et al.,¹⁷ who design a fully convolutional network that incorporates dilated convolution. This technique significantly widens the receptive field of the network, thus enabling it to capture more contextual information, which translates to improved image quality. Wu et al.¹⁹ attempted to learn the results of the WLS model in low-resolution (LR) space, then upsample the LR results to high-resolution with the trainable guided filtering layer. Thanks to the downsample-process-upsample strategy, the method significantly improves the computational efficiency. However, a common limitation among these learning-based approaches^18,17,19 is their lack of support for parameter tuning, which makes them less adaptable to various image processing tasks. To overcome this limitation, Yang et al.²⁰ proposed an unsupervised learning framework to approximate the results of the WLD model. Different from existing supervised learning-based methods, which require a huge amount of training data with ground truth labels, unsupervised training can be less costly. Nevertheless, the fixed structure of neural networks employed within this framework confines the solution space, which consequently affects the accuracy of the approximations.

In summary, while the original WLS solution delivers promising results, its computational demands are significant. The 1D approximations^15,16 of the WLS solution are computationally efficient but often introduce streak artifacts. Furthermore, existing deep learning-based approximations^18,17,19,20 show potential for the task, but they frequently suffer from limited flexibility or compromised output quality. Therefore, it is imperative to develop a fast solution to the WLS model that not only produces high-fidelity results but is also highly efficient.

Proposed method

Theoretical background

The WLS model¹⁴ is given as follows:

\underset{S}{\arg min} \sum_{i} ((S_{i} - P_{i})^{2} + λ \sum_{* \in {h, v}} w_{*, i} (\partial_{*} S_{i})^{2})

(1)

where P and S are the input and output images, respectively. i is the pixel index.

\partial_{*} S_{i}

represents the gradients of S at pixel i in the horizontal (h) and vertical (v) directions. λ is the balancing parameter.

w_{*, i}

stands for the weight at pixel i, and it is

w_{*, i} = ((\partial_{*} P_{i})^{α} + ϵ)^{- 1}

(2)

where

ϵ

is a small constant that avoids division by 0. The matrix notation of equation (1) can be written as follows:

\underset{S}{\arg min} ‖ S - P ‖_{2}^{2} + λ \sum_{* \in {h, v}} S^{T} D_{*}^{T} W_{*} D_{*} S

(3)

where

‖ \cdot ‖_{2}

is the

L_{2}

norm.

D_{*}

is the discrete differentiation matrix,

W_{*}

is the diagonal matrix containing the weights

w_{*, i}

. The superscript T is matrix transpose. The minimum of equation (3) is achieved when its derivative to S equals 0, that is,

(S - P) + λ \sum_{* \in {h, v}} D_{*}^{T} W_{*} D_{*} S = 0

(4)

Then, we have the naive solution as follows:

S = {(I + λ \sum_{* \in {h, v}} D_{*}^{T} W_{*} D_{*})}^{- 1} P

(5)

where I is the identity matrix.

\sum_{* \in {h, v}} D_{*}^{T} W_{*} D_{*}

is the inhomogeneous Laplacian matrix. Equation (5) is a large sparse linear system that can be time-consuming to solve.

Proposed solution

Different from existing fast solutions that either explore 1D approximation or deep learning-based approximation, in this article, we focus on the WLS model defined over the 2D space. We decouple the gradient regularization from the weighted scheme so that we would have a constrained optimization problem. Then we introduce the augmented Lagrangian multipliers to transform it back into an unconstrained but decoupled optimization problem, which can be optimized alternatively. Each iteration of our solution consists of a least squares problem, together with two elementwise updates. Due to the fact that the Laplacian matrix in the least squares problem is homogeneous, the problem can be efficiently solved in the Fourier domain.

Specifically, we first introduce an auxiliary variable $Z_{*} = D_{*} S$ to decouple the gradient regularization from the weighted scheme. Then we can construct a constrained optimization problem from equation (3) as follows:

\underset{S, Z_{*}}{\arg min} \begin{matrix} } V e r t S - P ‖_{2}^{2} + λ \sum_{* \in {h, v}} Z_{*}^{T} W_{*} Z_{*} \\ s.t. Z_{*} = D_{*} S \end{matrix}

(6)

We construct the augmented Lagrange function as follows:

\begin{aligned} \underset{S, Z_{*}, A_{*}}{\arg min} ‖ S - P ‖_{2}^{2} + λ \sum_{* \in {h, v}} Z_{*}^{T} W_{*} Z_{*} \\ + \sum_{* \in {h, v}} (\frac{η}{2} ‖ D_{*} S - Z_{*} ‖_{2}^{2} + A_{*}^{T} (D_{*} S - Z_{*})) \end{aligned}

(7)

Let

U_{*} = (1 / η) A_{*}

, we have the following equation:

\begin{aligned} \sum_{* \in {h, v}} (\frac{η}{2} ‖ D_{*} S - Z_{*} ‖_{2}^{2} + A_{*}^{T} (D_{*} S - Z_{*})) \\ = \sum_{* \in {h, v}} (\frac{η}{2} ‖ D_{*} S - Z_{*} + U_{*} ‖_{2}^{2} - \frac{η}{2} ‖ U_{*} ‖_{2}^{2}) \end{aligned}

(8)

Then equation (7) can be rewritten as follows:

\begin{aligned} \underset{S, Z_{*}, U_{*}}{\arg min} ‖ S - P ‖_{2}^{2} + λ \sum_{* \in {h, v}} Z_{*}^{T} W_{*} Z_{*} \\ + \sum_{* \in {h, v}} (\frac{η}{2} ‖ D_{*} S - Z_{*} + U_{*} ‖_{2}^{2} - \frac{η}{2} ‖ U_{*} ‖_{2}^{2}) \end{aligned}

(9)

Equation (9) can be solved by alternatively optimizing it with respect to S,

Z_{*}

, and

U_{*}

, then we have the following subproblems:

\underset{S}{\arg min} ‖ S - P ‖_{2}^{2} + \sum_{* \in {h, v}} \frac{η}{2} ‖ D_{*} S - Z_{*} + U_{*} ‖_{2}^{2}

(10)

\underset{Z_{*}}{\arg min} \sum_{* \in {h, v}} (λ Z_{*}^{T} W_{*} Z_{*} + \frac{η}{2} ‖ D_{*} S - Z_{*} + U_{*} ‖_{2}^{2})

(11)

\underset{U_{*}}{\arg min} \sum_{* \in {h, v}} (\frac{η}{2} ‖ D_{*} S - Z_{*} + U_{*} ‖_{2}^{2} - \frac{η}{2} ‖ U_{*} ‖_{2}^{2})

(12)

Equation (10) can be efficiently solved in the Fourier domain as follows:

S^{{k + 1}} = F^{- 1} (\frac{2 F (P) + η \sum_{* \in {h, v}} \bar{F (D_{*})} \circ F (Z_{*}^{{k}} - U_{*}^{{k}})}{2 F (1) + η \sum_{* \in {h, v}} \bar{F (D_{*})} \circ F (D_{*})})

(13)

where the superscript

{k}

indicates the k-th iteration.

F (x)

and

F^{- 1} (x)

represent the fast Fourier transform (FFT) and the inverse fast Fourier transform (IFFT) of x, respectively.

\bar{F (x)}

is the complex conjugate of

F (x)

F (1)

denotes the FFT of the delta function.

\circ

is the elementwise product. It should be emphasized that the division operation in equation (13) is also elementwise. Therefore, equation (13) can be evaluated efficiently.

Equation (11) can be conveniently solved as follows:

Z_{*}^{{k + 1}} = η (2 λ W_{*} + η I)^{- 1} (D_{*} S^{{k + 1}} + U_{*}^{{k}})

(14)

Equation (14) involves a matrix inverse operation, but it is efficient because that the matrix

2 λ W_{*} + η I

is a diagonal matrix. The inverse of the diagonal matrix can be simply calculated by inverting its diagonal elements. In fact, equation (14) can be implemented as elementwise operations, that is,

Z_{*}^{{k + 1}} = diag {(η (2 λ W_{*} + η I)}^{- 1}) \circ (D_{*} S^{{k + 1}} + U_{*}^{{k}})

, where

diag (x)

is the vector formed by the diagonal elements of x.

Equation (12) can be optimized as follows:

U_{*}^{{k + 1}} = U_{*}^{{k}} + D_{*} S^{{k + 1}} - Z_{*}^{{k + 1}}

(15)

Building upon the mathematical derivation, we sketch the algorithm of the proposed method as Algorithm 1. Our algorithm is iterative, where each iteration consists of three steps corresponding to equations (13) to (15). The main computational load in each iteration involves a few FFT and IFFT operations, which take O

(n \log n)

time, where n is the number of pixels. Considering there are K iterations, the total running time is O

(K n \log n)

. Note that for efficiency, we precompute the FFTs that remain constant during the iterative process. Therefore, in each iteration, we only need to calculate two FFTs and one IFFT.

Algorithm 1

The proposed method

Input: image P, parameters λ, α, iteration count K;

1: Initialize

S^{{0}} = P

Z_{*} = 0

U_{*} = 0

, where

* \in {h, v}

;

2: Precompute

F (P), F (\partial_{*}), \bar{F (\partial_{*})}

;

3: Precompute the denominator of equation (13);

4: for

k = 0 : K - 1

5: With

Z_{*}^{{k}}

and

U_{*}^{{k}}

, solve equation (13) for

S^{{k + 1}}

;

6: With

S^{{k + 1}}

and

U_{*}^{{k}}

, solve equation (14) for

Z^{{k + 1}}

;

7: With

S^{{k + 1}}

and

Z_{*}^{{k + 1}}

, solve equation (15) for

U^{{k + 1}}

;

8: end for

Output: image

S^{{K}}

;

Akin to the original WLS model, there are three tunable parameters in our method, that is, λ, α, and K. The former two control the smoothness of the resultant images, while the third one is used to control the approximation accuracy. The recommended ranges for these parameters are $λ \in [0.001, 0.1]$ , $α \in [1, 4]$ , $K = 10$ . Empirically, these values would produce satisfactory results.

Applications

Due to the limitations of existing medical imaging systems, medical images may be subjected to poor visibility and insufficient details, thus necessitating detail and tone enhancement. The proposed method can facilitate both of the two tasks. In this section, we will detail the procedure of applying the proposed method to these two tasks.

Detail enhancement

For convenient and clear descriptions, we denote the function provided by Algorithm 1 as $f (\cdot)$ . The procedure of detail enhancement starts by deriving the base layer $S_{d}$ using the proposed method, that is,

S_{d} = f (P_{d})

(16)

where

P_{d}

is the input. Then, the detail layer can be derived by subtracting the base layer from the input image, that is,

δ_{d} = P_{d} - S_{d}

(17)

Ultimately, the detail layer is magnified m times and added back to the input to derive the final detail-enhanced image

O_{d}

, that is,

O_{d} = m δ_{d} + P_{d}

(18)

The detail-enhancing process is also illustrated in Figure 1. An example of a chest X-ray image detail enhancement process can be found in Figure 2. We highlight the mathematical symbols in blue to facilitate connections between the figures and the equations.

Figure 1.

The pipeline of the detail enhancement scheme.

Figure 2.

An example of the detail enhancement process: (a) input, (b) base layer, (c) detail layer, and (d) detail-enhanced image.

Figure 3 demonstrates a visual comparison between the proposed method and popular detail enhancement techniques, including histogram equalization and bilateral filtering. It is evident that histogram equalization fails to significantly enhance the details of the image. The bilateral filter enhances the details, but it suffers from the halo artifacts (see the alphabet A). In contrast, our method properly enhances the details while being free of the artifacts.

Figure 3.

Visual comparison to popular detail enhancement methods on (a) an input X-ray image. The details in the result of (b) histogram equalization are not significantly enhanced. The result of (c) bilateral filter suffers from the halo artifacts (see the alphabet A). In contrast, the details are properly enhanced in (d) our result without suffering from artifacts.

Tone enhancement

According to the retinex theory,²⁵ the observed image $P_{t}$ can be decomposed into the reflectance $O_{t}$ and illumination $L_{t}$ , that is,

P_{t} = O_{t} \circ L_{t}

(19)

The objective of tone enhancement is to estimate

O_{t}

from

P_{t}

. However, in order to do that, we need to first estimate the illumination map

L_{t}

. A coarse guess for

L_{t}

is the maximum intensity in all the three channels of

P_{t}

, that is,

\hat{L_{t}} = max_{⋆ \in {R, G, B}} P_{t}^{⋆} + ϵ

(20)

where

ϵ

is a compensation parameter and it is set as 0.2. The illumination map is refined with our method, which is then followed by gamma correction, that is,

L_{t} = f (\hat{L_{t}})^{γ}

(21)

Afterwards, we can estimate

O_{t}

from equation (19). It should be noted that

O_{t}

may contain noise and thus needs to be denoised and recomposed back with the input

P_{t}

O_{t}^{†} = (1 - L_{t}) \circ τ (O_{t}) + P_{t}

(22)

The tone enhancement process is also shown in Figure 4. An example of a retinal fundus image tone enhancement process can be found in Figure 5.

Figure 4.

The pipeline of the tone enhancement method.

Figure 5.

An example of the tone enhancement process: (a) input image, (b) initial illumination map, (c) refined illumination map, and (d) tone-enhanced image.

The visual comparison on the task of tone enhancement is shown in Figure 6. It can be observed that histogram equalization struggles with insufficient contrast. While the bilateral filter addresses the contrast issue, it lacks fine details. Conversely, our method produces a vivid tone-enhanced result with adequate contrast and details.

Figure 6.

Visual comparison to popular tone enhancement methods on (a) an input fundus image. The result of (b) histogram equalization lacks contrast. The result of (c) bilateral filter lacks details. Conversely, (d) our result is of high quality.

Experiment settings

We have experimented with the proposed method on the detail and tone enhancement of medical images. For the former one, we use the COVID-19 chest X-ray image dataset²⁶ (317 X-ray images) and a subset (243 CT and MRI images) from the CHAOS dataset.²⁷ For the latter one, we leverage the retinal fundus image dataset named CHASE-DB1,²⁸ which consists of 28 retinal images collected from 14 children. The proposed method is rigorously evaluated using the test sets from these datasets.

To comprehensively evaluate the proposed method, we compare it to a variety of fast WLS methods, including the WLS,¹⁴ fast global WLS (FG-WLS),¹⁶ semi-global WLS (SG-WLS),¹⁵ deep edge-aware WLS (DE-WLS),¹⁸ context aggregation WLS (CA-WLS)¹⁷ deep guided WLS (DG-WLS),¹⁹ and unsupervised WLS (UN-WLS).²⁰ It should be noted that for traditional methods, the parameters are best tuned, whereas for deep learning-based methods, we use the officially pretrained model, where the parameters are not tunable.

We employ both the perception-based and distortion-based image quality measures to quantitatively assess the proposed method. The perception-based measures for the task of detail enhancement include the spatial-spectral entropy-based quality (SSEQ)²⁹ and convolutional neural networks image quality assessment (CIQA).³⁰ The perception-based measures for tone enhancement are tone mapped image quality index (TMQI),³¹ along with its subindices, that is, fidelity and naturalness. The distortion-based measures for both tasks are the peak signal-to-noise ratio (PSNR) and structural similarity index (SSIM). Furthermore, we demonstrate visual results for qualitative evaluation of the two tasks. It should be noted that for SSEQ and CIQA, smaller values indicate better quality, while for TMQI, PSNR, and SSIM, larger values correspond to better quality. Finally, we record and compare the running times of various methods.

Results

Detail enhancement

The perception-based quality metrics of detail enhancement on the X-ray dataset and CHAOS subset are shown in Table 1 and Figure 7, respectively. It can be drawn that our method consistently achieves a better perceptual quality (including both SSEQ and CIQA) than existing fast WLS methods across all the datasets involved. In fact, our method delivers a comparable perceptual quality to that of the original WLS method. This can be explained by the distortion-based quality metrics of detail enhancement shown in Figure 8, where we see that for different sets of parameters, our method could accurately approximate the results of the original WLS method. Specifically, considering there are 10 iterations (i.e. the value that we use throughout the article), our method could achieve PSNR $\geq 39.39$ , while SSIM $\geq 0.99$ , which is sufficiently accurate for real applications. Furthermore, we observe that, as the number of iterations increases, the approximation accuracy gets even better.

Figure 7.

The perception-based quality measures of detail enhancement on the CHAOS subset²⁷ (CT and MRI images). Smaller values of SSEQ and CIQA indicate better quality. CT: computer tomography; MRI: magnetic resonance images; SSEQ: spatial-spectral entropy-based quality; CIQA: convolutional neural networks image quality assessment.

Figure 8.

The distortion-based quality measures for detail enhancement on the X-ray dataset.²⁶ Larger values of PSNR and SSIM indicate better quality. PSNR: peak signal-to-noise ratio; SSIM: structural similarity index.

Table 1.

The perception-based quality measures of detail enhancement on the COVID-19 X-ray dataset.²⁶

Method	Covid		Normal		Viral
Method	SSEQ	CIQA	SSEQ	CIQA	SSEQ	CIQA
WLS	21.16	17.59	29.49	25.20	31.36	25.21
FG-WLS	25.35	25.77	33.58	41.04	34.90	41.43
SG-WLS	25.25	25.41	34.53	40.26	35.02	41.42
DE-WLS	25.79	24.70	35.62	38.83	35.94	40.49
CA-WLS	24.71	19.68	31.86	32.66	31.53	31.81
DG-WLS	34.28	20.88	32.22	27.47	32.54	28.74
UN-WLS	24.18	25.07	35.02	40.21	35.59	41.34
Ours	21.15	17.67	29.46	25.60	31.36	26.01

Bold: best, Underline: second best. Smaller values of SSEQ and CIQA indicate better quality.

SSEQ: spatial-spectral entropy-based quality; CIQA: convolutional neural networks image quality assessment; WLS: weighted least squares; FG-WLS: fast global weighted least squares; SG-WLS: semi-global weighted least squares; DE-WLS: deep edge-aware weighted least squares; CA-WLS: context aggregation weighted least squares; DG-WLS: deep guided weighted least squares; UN-WLS: unsupervised weighted least squares.

The superiority of our method on detail enhancement can be also validated from the qualitative results on an X-ray image shown in Figure 9. We can see that due to the separated 1D implementation, FG-WLS suffers from severe streak artifacts (check the areas above the pelvis). UN-WLS does not suffer from artifacts, but its result is of insufficient details. Furthermore, it tends to deviate the color of the image. On the other hand, our method is able to derive a promising result, which is comparable to that of the original WLS method. In the meantime, our method can be much more efficient than WLS. This advantage is further supported by the qualitative results on a CT image presented in Figure 10.

Figure 9.

Qualitative results of detail enhancement on (a) an input X-ray image. The result of (b) FG-WLS suffers from the streak artifacts (see the areas above the pelvis). The result of (c) UN-WLS contains insufficient details. (e) Our result is as good as that of (d) WLS, but it can be derived in a much more efficient way. WLS: weighted least squares; FG-WLS: fast global weighted least squares; UN-WLS: unsupervised weighted least squares.

Figure 10.

Qualitative results of detail enhancement on (a) an input CT image. The result of (b) FG-WLS suffers from the gradient reversal artifacts (around salient edges). The result of (c) UN-WLS lacks details. (e) Our result is of high quality akin to the result of (d) WLS, but is achieved in a significantly more efficient manner. CT: computer tomography; WLS: weighted least squares; FG-WLS: fast global weighted least squares; UN-WLS: unsupervised weighted least squares.

Tone enhancement

The perception-based quality metrics of tone enhancement are shown in Table 2. It can be observed that our method outperforms existing fast WLS methods in terms of perceptual quality, including TMQI and its subindices. As a matter of fact, our method delivers a comparable perceptual quality to that of the original WLS method. This can be observed from the distortion-based quality metrics illustrated in Figure 11, where our method accurately approximates the results of the original WLS method for various sets of parameters. Specifically, considering there are 10 iterations (i.e. the value that we use throughout the article), our method achieves PSNR $\geq 50.09$ , while SSIM $\geq 0.99$ , which is deemed sufficiently accurate for practical applications. Moreover, the accuracy gets even higher as we increase the number of iterations.

Figure 11.

The distortion-based quality measures for tone enhancement. Larger values of PSNR and SSIM indicate better quality. PSNR: peak signal-to-noise ratio; SSIM: structural similarity index.

Table 2.

The perception-based quality measures on tone enhancement.

Method	TMQI	Fidelity	Naturalness
WLS	0.7911	0.9329	0.0086
FG-WLS	0.7779	0.8907	0.0049
SG-WLS	0.7827	0.9049	0.0066
DE-WLS	0.7872	0. 9186	0.0082
CA-WLS	0.7861	0.9158	0.0077
DG-WLS	0.7848	0.9123	0.0070
UN-WLS	0.7864	0.9136	0.0093
Ours	0.7905	0.9300	0.0090

Bold: best, Underline: second best. Larger values of TMQI, fidelity, and naturalness indicate better quality.

TMQI: tone mapped image quality index; WLS: weighted least squares; FG-WLS: fast global weighted least squares; SG-WLS: semi-global weighted least squares; DE-WLS: deep edge-aware weighted least squares; CA-WLS: context aggregation weighted least squares; DG-WLS: deep guided weighted least squares; UN-WLS: unsupervised weighted least squares.

The qualitative results shown in Figure 12 further corroborate the advantages of our method. It is evident that FG-WLS suffers from noticeable streak artifacts (as seen in the top-right corner of the circle). While UN-WLS does not exhibit artifacts, its output lacks contrast and details. Conversely, our method yields an impressive result, which is almost identical to the original WLS method. At the same time, our method is notably more efficient than WLS.

Figure 12.

Qualitative results of tone enhancement on (a) an input fundus image. The result of (b) FG-WLS is affected by the streak artifacts (see the top-right corner of the circle). The result of (c) UN-WLS lacks contrast and details. (e) Our result is equally impressive as that of (d) WLS while being much faster to derive. WLS: weighted least squares; FG-WLS: fast global weighted least squares; UN-WLS: unsupervised weighted least squares.

Efficiency

The objective of this article is to develop an efficient solution to the WLS model. Therefore, the running time is a crucial evaluation metric in this article. To this end, we implement our method on both central processing unit (CPU) and graphics processing unit (GPU), and compare its running time to that of existing methods. The running times reported are recorded on an AMD Ryzen 4800H, an NVIDIA RTX 3070 (8 GB GDDR6), and 16 GB main memory. Note that the time is measured in seconds. The running time comparison is shown in Table 3. We can see that our method is faster than most of the existing WLS solutions, except DG-WLS. However, the high accuracy of DG-WLS is attributed to the downsampling-processing-upsampling strategy, which processes the image at LR, and then upsamples it with the assistance of a high-resolution input image. The strategy would significantly accelerate image processing, but it would also sacrifice the image quality, as demonstrated in the previous sections. Our method is highly efficient, it is able to process 720P color images at more than 13 frames per second. Compared to the original WLS, our CPU implementation is about 10 to 14 times faster, our GPU implementation is about 28 to 142 times faster. This validates the fact that our method is a fast solution to the WLS model.

Table 3.

Running times of various methods on images with different resolutions (measured in seconds).

Method	$320 \times 240$	$640 \times 480$	$1280 \times 720$
WLS	0.6388	3.0029	10.3171
FG-WLS	0.0895	0.3775	1.2411
SG-WLS	0.1259	0.5103	1.6868
DE-WLS	0.1815	0.5994	1.7010
CA-WLS	0.0080	0.0310	0.0870
DG-WLS	0.0030	0.0100	0.0290
UN-WLS	0.0058	0.0295	0.0770
Ours (CPU)	0.0454	0.2971	0.8114
Ours (GPU)	0.0257	0.0453	0.0725

WLS: weighted least squares; FG-WLS: fast global weighted least squares; SG-WLS: semi-global weighted least squares; DE-WLS: deep edge-aware weighted least squares; CA-WLS: context aggregation weighted least squares; DG-WLS: deep guided weighted least squares; UN-WLS: unsupervised weighted least squares; CPU: central processing unit; GPU: graphics processing unit.

Parameters

The parameters λ, α, and K have a significant impact on the enhanced results, as demonstrated in Figure 13. Specifically, the balancing parameter λ and the edge amplitude α control the spatial scale of the elements to be enhanced, that is, larger values of λ and α render more elements to be recognized as details, thus enabling the proposed method to enhance larger-scale details. The iteration count K determines the accuracy in approximating the original WLS solution. A larger value of K corresponds to a better approximation accuracy. Empirically, the enhanced images would barely change when $K > 10$ .

Figure 13.

The impact of parameters on the enhanced results. λ: balancing parameter; α: edge amplitude; K: iteration count.

Discussion

Our method achieves impressive results as demonstrated in the previous sections. However, the performance can be further improved by tuning the best parameters for individual images. An important common task for various applications is automatic parameter setting. To this end, the images with the best quality need to be chosen through meticulously designed subjective and objective studies. This can be done by gathering a diverse set of images and having them rated by professionals to ensure a comprehensive understanding of what constitutes “best” in different applications. Then the relationship between the optimal parameters and the image features, such as brightness, contrast, texture, and color distribution can be established with supervised learning. The importance of features can be determined according to the perception-based image quality. By understanding how these features correlate with the optimal parameters, we can develop models or algorithms that predict the parameter settings for new images based on their features.

In this article, we experiment with the proposed method on the detail enhancement of the X-ray images and the tone enhancement of the retinal fundus image. For the X-ray images, our focus is on enhancing fine details that are crucial for accurate diagnosis. For the retinal fundus images, which capture the inner surface of the eye, including the retina, optic disc, macula, and blood vessels, our method focuses on improving the tone and contrast. These images are essential for diagnosing and monitoring various ocular diseases, such as diabetic retinopathy, glaucoma, and age-related macular degeneration. However, our method is not limited to only detail enhancement or tone enhancement. As a matter of fact, these two enhancement techniques can be used in combination to improve the overall quality of the images, thereby providing clearer and more detailed visual information for medical practitioners.

Enhanced image quality can lead to better visualization of critical structures and pathologies, which results in more accurate diagnoses and reduces chances of misdiagnosis. Furthermore, improved diagnostic capabilities can lead to earlier detection and treatment of conditions, ultimately enhancing patient outcomes and satisfaction. The proposed method can be integrated into clinical workflows by developing plugins or modules with user-friendly interfaces and embedding them into radiology software. The enhanced images can be also linked to the electronic health records, which allows for seamless access to improved images alongside patient history and other relevant data. However, technical challenges may arise when integrating the new method with existing imaging systems and workflows as it requires significant IT support and resources. Moreover, the proposed method must undergo rigorous validation to ensure its safety and efficacy in clinical settings. Obtaining regulatory approval can be a lengthy and complex process.

Limitations

While the proposed method for detail and tone enhancement in medical images shows promising results, several limitations need to be taken into account:

The effectiveness of the proposed method may diminish when applied to images with extreme noise levels. Medical images may contain significant noise due to various factors such as low signal-to-noise ratios, patient movement, or equipment limitations. In these situations, the enhancement process may inevitably amplify the noise rather than the relevant features, thus leading to degraded image quality. This limitation suggests that additional preprocessing steps, such as advanced denoising techniques, may be necessary to improve the performance of the proposed method on images with heavy noise.

The performance of the proposed method is likely to be compromised when dealing with very LR images. LR images inherently lack detail, which can impede the enhancement process. The proposed method relies on the presence of sufficient image information to effectively enhance details and tones. In scenarios where the resolution is too low, the method may struggle to produce meaningful enhancements, thus resulting in images that appear blurred or lack clarity.

Although the method has been validated on specific types of medical images, such as COVID-19 chest X-rays, MRI, CT, and retinal fundus images, its applicability to other imaging modalities (e.g. ultrasound and endoscopy images) or non-medical images remains uncertain. Different imaging techniques possess unique characteristics, such as varying noise profiles, contrast levels, and types of artifacts. As a result, the performance of the proposed method may not be consistent across all image types.

Conclusion and future work

In this article, we propose an efficient WLS method for the detail and tone enhancement of the medical images. To this end, we propose to firstly decouple the gradient regularization from the weighted scheme so that the foundation of the inhomogeneous Laplacian matrix is broken. Then we propose a novel solution based on the ADMM and the Fourier domain optimization. Different from existing deep learning-based approximations of the WLS method, our method support parameter tuning, which are more flexible. Different from existing separated 1D implementations, our method focuses on 2D image processing, thus being free of the streak artifacts. The proposed method is highly efficient, it is able to process 720P color images in real time. To evaluate the proposed method, we experiment with the proposed method on the detail enhancement of the X-ray images and the tone enhancement of the retinal fundus images. A variety of perception- and distortion-based image quality measures are leveraged, including SSEQ, CIQA, TMQI, PSNR, and SSIM. Both quantitative and qualitative results substantiate the superiority of the proposed method.

There are several promising directions for future research and applications. The enhanced image quality produced by our method could be integrated into automated diagnostic systems to improve the accuracy of image analysis and interpretation. For instance, by providing clearer and more detailed images, the method could assist radiologists and ophthalmologists in diagnosing conditions such as diabetic retinopathy, glaucoma, and other ocular diseases more effectively. Furthermore, the versatility of the proposed method allows for its potential application across various medical imaging modalities, including endoscopy, and ultrasound images. Future research endeavors will focus on extending the applicability of our method to additional imaging modalities. Besides, creating user-friendly software tools that incorporate our method could facilitate its adoption by medical professionals. Such tools could provide intuitive interfaces for image enhancement, which allows clinicians to easily improve image quality without requiring extensive technical knowledge. Finally, even though the proposed solution is highly efficient, it is not yet ready for extremely high-resolution images. Further efforts are needed to improve the solution. Possible directions include the moth-flame,^32,33 cuckoo,^34,35 multi-verse,^36,37 and chimp^38,39 optimization algorithms.

Footnotes

Acknowledgements

We would like to express our deepest gratitude to the people who developed the dataset involved in this article, including Pranav Raikote and Christopher G. Owen, etc.

Contributorship

YY and WB: conceptualization; WD, WB, and YY: data curation; WB, WD, and YY: methodology; WB, QZ, QG, SB, and YY: formal analysis; YY: supervision; WB, QZ, QG, SB, and YY: validation; WB: writing–original draft; WD, QZ, QG, SB, and YY: writing–review and editing.

Consent statement

Patient consent was not required for the present manuscript, as the datasets in this study are publicly available.

Declaration of conflicting interests

The author(s) declared no potential conflicts of interest with respect to the research, authorship, and/or publication of this article.

Ethics approval

Not applicable.

Funding

The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: The author(s) disclosed receipt of the following financial support for the research, authorship, and/or publication of this article: This project was funded by National Natural Science Foundation of China (Ref: 61402205), China Postdoctoral Science Foundation (Ref: 2015M571688), Jiangsu University (Ref: 13JDG085), Jiangsu Province Elderly Health Project (Ref: LD2021039), and Zhenjiang Key Laboratory of Health and Life Sciences (Ref: GZSYS202304).

Guarantor

YY.

ORCID iD

Yang Yang

References

Doi

. Diagnostic imaging over the last 50 years: research and development in medical imaging science and technology. Phys Med Biol 2006; 51: R5.

Philip

Cowie

Olson

. The impact of the Health Technology Board for Scotland’s grading model on referrals to ophthalmology services. Br J Ophthalmol 2005; 89: 891–896.

Macgillivray

Cameron

Zhang

, et al. Suitability of UK Biobank retinal images for automatic analysis of morphometric properties of the vasculature. PLoS ONE 2015; 10: e0127914.

Kim

. Contrast enhancement using brightness preserving bi-histogram equalization. IEEE Trans Consum Electron 1997; 43: 1–8.

Roy

Bhalla

Patel

. Mathematical analysis of histogram equalization techniques for medical image enhancement: a tutorial from the perspective of data loss. Multim Tools Appl 2024; 83: 14363–14392.

Abdullah-Al-Wadud

Kabir

Dewan

MAA

, et al. A dynamic histogram equalization for image contrast enhancement. IEEE Trans Consumer Electron 2007; 53: 593–600.

Veluchamy

Bhandari

Subramani

. Optimized Bezier curve based intensity mapping scheme for low light image enhancement. IEEE Trans Emerg Top Comput Intell 2022; 6: 602–612.

Awarayi

Twum

Hayfron-Acquah

, et al. A bilateral filtering-based image enhancement for Alzheimer disease classification using CNN. PLoS ONE 2024; 19: 1–15.

Rajendran

Rao

Agaian

, et al. A versatile edge preserving image enhancement approach for medical images using guided filter. In: IEEE international conference on systems, man, and cybernetics, Budapest, Hungary, October 9-12, 2016, pp. 2341–2346. IEEE.

10.

Bian

Yang

. Fast bilateral weighted least square for the detail enhancement of COVID-19 chest X-rays. Digit Health 2023; 9: 20552076231200981.

11.

Chen

Pawar

Ekanayake

, et al. Deep learning for image enhancement and correction in magnetic resonance imaging—state-of-the-art and challenges. J Digit Imaging 2023; 36: 204–230.

12.

Alla Takam

Samba

Tchagna Kouanou

, et al. Spark architecture for deep learning-based dose optimization in medical imaging. Inform Med Unlocked 2020; 19: 100335.

13.

Ghandour

El-Shafai

El-Rabaie

. Medical image enhancement algorithms using deep learning-based convolutional neural network. J Opt 2023; 52: 1931–1941.

14.

Farbman

Fattal

Lischinski

, et al. Edge-preserving decompositions for multi-scale tone and detail manipulation. ACM Trans Graph 2008; 27: 67.

15.

Liu

Chen

Shen

, et al. Semi-global weighted least squares in image filtering. In: IEEE international conference on computer vision, Venice, Italy, October 22-29, 2017, pp. 5862–5870. IEEE Computer Society.

16.

Min

Choi

, et al. Fast global image smoothing based on weighted least squares. IEEE Trans Image Process 2014; 23: 5638–5653.

17.

Chen

Koltun

. Fast image processing with fully-convolutional networks. In: IEEE international conference on computer vision, Venice, Italy, October 22-29, 2017, pp.2516–2525. IEEE Computer Society.

18.

Ren

JSJ

Yan

, et al. Deep edge-aware filters. In: Proceedings of the 32nd international conference on machine learning, Vol. 37, 2015, pp.1669–1678. JMLR.org

19.

Zheng

Zhang

, et al. Fast end-to-end trainable guided filter. In: IEEE conference on computer vision and pattern recognition, Salt Lake City, UT, USA, June 18-22, 2018, pp. 1838–1847. Computer Vision Foundation / IEEE Computer Society.

20.

Yang

Zeng

, et al. Weighted least square filter via deep unsupervised learning. Multim Tools Appl 2024; 83: 31361–31377.

21.

Alenezi

Santosh

. Geometric regularized hopfield neural network for medical image enhancement. Int J Biomed Imaging 2021; 2021: 6664569:1.

22.

Liu

, et al. Structure and illumination constrained GAN for medical image enhancement. IEEE Trans Med Imaging 2021; 40: 3955–3967.

23.

Zhu

You

, et al. Pre-trained diffusion models for plug-and-play medical image enhancement. In: Medical Image Computing and Computer Assisted Intervention - MICCAI 2023 - 26th International Conference, Vancouver, BC, Canada, October 8-12, 2023, Proceedings, Part II, series Lecture Notes in Computer Science, 14222, pp. 3–13. Springer.

24.

Perona

Malik

. Scale-space and edge detection using anisotropic diffusion. IEEE Trans Pattern Anal Mach Intell 1990; 12: 629–639.

25.

Guo

Ling

. LIME: low-light image enhancement via illumination map estimation. IEEE Trans Image Process 2017; 26: 982–993.

26.

Raikote

COVID-19 image dataset

https://www.kaggle.com/datasets/pranavraikokte/covid19-image-dataset , 2024.

27.

Kavur

Selver

Dicle

, et al. CHAOS Challenge - combined CT-MR healthy abdominal organ segmentation. Medical Image Anal 2021; 69: 101950.

28.

Owen

Rudnicka

Mullen

, et al. Measuring retinal vessel tortuosity in 10-year-old children: validation of the computer-assisted image analysis of the retina (CAIAR) program. Invest Ophthalmol Vis Sci 2009; 50: 2004–2010.

29.

Liu

Huang

, et al. No-reference image quality assessment based on spatial and spectral entropies. Signal Process Image Commun 2014; 29: 856–863.

30.

Kang

, et al. Convolutional neural networks for no-reference image quality assessment. In: IEEE conference on computer vision and pattern recognition, CVPR 2014, Columbus, OH, USA, June 23-28, 2014, pp. 1733–1740. IEEE Computer Society.

31.

Yeganeh

Wang

. Objective quality assessment of tone-mapped images. IEEE Trans Image Process 2013; 22: 657–667.

32.

SMoth-flame optimization algorithm Mirjalili. A novel nature-inspired heuristic paradigm. Knowl Based Syst 2015; 89: 228–249.

33.

Shehab

Abualigah

Hamad

, et al. Moth-flame optimization algorithm: variants and applications. Neural Comput Appl 2020; 32: 9859–9884.

34.

Yang

Deb

. Cuckoo search via Lévy flights. In: World congress on nature & biologically inspired computing, NaBIC 2009, 9-11 December 2009, Coimbatore, India, pp. 210–214}. IEEE.

35.

Shehab

Khader

Laouchedi

, et al. Hybridizing cuckoo search algorithm with bat algorithm for global numerical optimization. J Supercomput 2019; 75: 2395–2422.

36.

Mirjalili

Hatamlou

. Multi-verse optimizer: a nature-inspired algorithm for global optimization. Neural Comput Appl 2016; 27: 495–513.

37.

Shehab

Abualigah

. Opposition-based learning multi-verse optimizer with disruption operator for optimization problems. Soft Comput 2022; 26: 11669–11693.

38.

Khishe

Mosavi

. Chimp optimization algorithm. Expert Syst Appl 2020; 149: 113338.

39.

Daoud

Shehab

Abualigah

, et al. Hybrid modified chimp optimization algorithm and reinforcement learning for global numeric optimization. J Bionic Eng 2023; 20: 2896–2915.