Method and system for identifying target lesion based on medical image

Abstract

The application relates to the technical field of medical image processing, and provides a target lesion identification method and system based on medical images, which solves the problems of limited and unstable lesion identification precision. The method comprises the following steps: acquiring a cone beam CT image of a mandible before oral implantation and reference gray data; performing metal artifact reduction and adaptive histogram equalization algorithm enhancement processing with limited contrast to obtain an enhanced image; extracting a three-dimensional gray histogram of an implant target area to calculate gray skewness and gray kurtosis, and extracting a texture feature set; performing principal component analysis on the texture feature set to obtain principal component features, and splicing the principal component features with the gray skewness and the gray kurtosis into a fusion feature vector; analyzing the vector by using a density peak value clustering algorithm, and mapping the result to an anatomical partition, and performing Z-score testing in combination with the reference gray data to identify a bone lesion area. The application can realize accurate and automatic positioning of a lesion.

Description

Technical Field

[0001] This application relates to the technical field of medical image processing, and in particular to a method and system for identifying target lesions based on medical images. Background Technology

[0002] In the preoperative planning of oral implant surgery, accurate assessment of the bone condition of the mandible and potential lesions is crucial. Furthermore, lesion identification based on medical imaging can help doctors objectively judge the implant conditions, which has a clear application prospect in improving the success rate of surgery and realizing personalized precision medicine.

[0003] Currently, clinical practice mainly relies on doctors' visual assessment and experience in cone-beam CT images. In addition, single grayscale threshold segmentation or basic texture statistical features are often used for region division and anomaly detection.

[0004] However, these existing methods have the following drawbacks: on the one hand, visual assessment is highly subjective and it is difficult to achieve stable and repeatable quantitative analysis; on the other hand, methods that rely on only a single feature or simple rules have limited ability to distinguish complex and variable trabecular bone structures and image artifact interference, ultimately resulting in insufficient recognition accuracy, low specificity, and low degree of automation. Summary of the Invention

[0005] The purpose of this application is to provide a method and system for identifying target lesions based on medical images, so as to solve the problems of limited accuracy and insufficient stability of the identification results of mandibular lesions in the prior art.

[0006] To address the aforementioned technical problems, in a first aspect, this application provides a method for identifying target lesions based on medical images, comprising:

[0007] Obtain cone-beam CT images of the mandible before oral implant surgery and reference grayscale data corresponding to the implantation target area;

[0008] The cone-beam CT images are processed using a metal artifact reduction algorithm to obtain artifact-suppressed images;

[0009] Based on the artifact-suppressed image, a contrast-limited adaptive histogram equalization algorithm is used to perform local contrast enhancement processing to obtain an enhanced image.

[0010] Extract the three-dimensional grayscale histogram of the implantation target region from the enhanced image, and calculate the grayscale skewness and grayscale kurtosis of the three-dimensional grayscale histogram;

[0011] Feature extraction is performed on the implantation target region in the enhanced image to obtain a texture feature set, wherein the texture feature set includes contrast, correlation, energy and homogeneity;

[0012] Principal component analysis is performed on the texture feature set to obtain principal component features that characterize the main variation direction of the trabecular bone structure. The gray-level skewness, gray-level kurtosis and the principal component features are concatenated according to preset weights to obtain a fused feature vector.

[0013] The fused feature vector is analyzed using a density peak clustering algorithm to obtain clustering results. The clustering results are then mapped to the anatomical partitions corresponding to the implantation target area. The Z-score test is performed on each sub-region corresponding to the clustering results in conjunction with the reference grayscale data to obtain the bone lesion area.

[0014] Optionally, the density peak clustering algorithm is used to analyze the fused feature vector to obtain clustering results. The clustering results are then mapped to the anatomical partitions corresponding to the implantation target region. A Z-score test is performed on each sub-region corresponding to the clustering results using the reference grayscale data to obtain the bone lesion region, including:

[0015] Based on the fused feature vector, the density peak clustering algorithm is used to calculate the local density value and relative distance value of each data point;

[0016] By combining the local density value and the relative distance value, a decision map is constructed. Based on the decision map, multiple candidate cluster centers are identified, and the clustering result is obtained according to the affiliation relationship between each data point and the candidate cluster centers.

[0017] Based on the clustering results, the implantation target region is divided into multiple clustered sub-regions, and each clustered sub-region is mapped to a corresponding anatomical partition. The Z-score is calculated by combining the mean gray value of the clustered sub-region with the mean gray value of the reference gray data.

[0018] Based on the distribution of the local density values ​​and the relative distance values ​​corresponding to the clustered sub-regions, a consistency metric is calculated within each clustered sub-region.

[0019] Combining the preset anatomical constraint rules corresponding to the anatomical partitions with the consistency metric, and introducing the structural tensor to optimize the boundaries of the corresponding clustered sub-regions, an optimized sub-region is generated. The preset anatomical constraint rules are based on the spatial location and morphology of the mandibular anatomical diagram, and the structural tensor is obtained based on image gradient calculation.

[0020] Calculate the distribution deviation between the grayscale distribution of the optimized sub-region and the normal grayscale model in the reference grayscale data;

[0021] The Z-score, the consistency metric, and the corresponding distribution deviation are weighted and summed to obtain the anomaly confidence of each optimized sub-region.

[0022] Candidate lesion regions with an abnormal confidence level exceeding a preset abnormal threshold are selected, and spatial aggregation is performed on the candidate lesion regions to obtain bone lesion regions.

[0023] Optionally, the step of weighted summing of the Z-score, the consistency metric, and the corresponding distribution deviation to obtain the anomaly confidence of each optimized sub-region includes:

[0024] The Z-score, the consistency metric, and the distribution deviation are input into a pre-trained three-layer fully connected neural network, which then predicts the fusion coefficients of each item.

[0025] Based on the fusion coefficient, the Z-score, the consistency metric, and the distribution deviation are weighted and summed to obtain a comprehensive anomaly index.

[0026] The composite anomaly index is converted into an intermediate confidence level using the Sigmoid function;

[0027] Based on the prior Gaussian distribution of normal gray values ​​in the reference gray data, the Bayesian formula is used to calculate the posterior probability of the intermediate confidence level to obtain the abnormal confidence level of each optimized sub-region.

[0028] Optionally, the step of processing the cone-beam CT image using a metal artifact reduction algorithm to obtain an artifact-suppressed image includes:

[0029] The original projection data of the cone-beam CT image is read to obtain the original sine wave. The sine wave is then segmented into metal regions to obtain a metal region mask.

[0030] Based on the metal region mask, a metal artifact reduction algorithm is used to extract multi-scale features from the original sine wave to obtain contextual features;

[0031] The feature sequence corresponding to the context features and the projection statistical features corresponding to the original sine curve are input into a gated recurrent unit for fusion to obtain an enhanced feature vector.

[0032] Based on the enhanced feature vector, data reconstruction is performed using a conditional generative adversarial network to obtain the repair data corresponding to the metal region mask.

[0033] The repaired data is fused with the original data corresponding to the original sine curve to obtain an intermediate sine curve. The intermediate sine curve is then optimized using a gradient descent iterative algorithm based on total variation minimization constraints to obtain the target sine curve.

[0034] The target sine wave is filtered and back-projected to reconstruct the image, resulting in an artifact-suppressed image.

[0035] Optionally, the step of performing local contrast enhancement processing based on the artifact-suppressed image using a contrast-limited adaptive histogram equalization algorithm to obtain an enhanced image includes:

[0036] Regional characteristic analysis was performed on the artifact-suppressed image to obtain a first characteristic region representing a high bone density region and a second characteristic region representing a low bone density region.

[0037] Based on the characteristic differences between the first characteristic region and the second characteristic region, regional statistical features are generated to guide the global enhancement strategy.

[0038] Based on the statistical characteristics of the region, a first constraint strategy, which mainly increases the clipping limit, is applied to the first characteristic region through a contrast-limited adaptive histogram equalization algorithm, and a second constraint strategy, which mainly decreases the clipping limit, is applied to the second characteristic region.

[0039] Meanwhile, based on the preset anatomical reference images corresponding to the first and second characteristic regions respectively, the execution process of the first and second constraint strategies is guided by a regularization term based on the gradient magnitude to obtain intermediate images.

[0040] The intermediate image is then subjected to boundary smoothing processing to obtain an enhanced image.

[0041] Secondly, this application provides a target lesion identification system based on medical images, comprising:

[0042] The acquisition module is used to acquire cone-beam CT images of the mandible before oral implant surgery and reference grayscale data corresponding to the implantation target area;

[0043] The processing module is used to process the cone-beam CT image using a metal artifact reduction algorithm to obtain an artifact-suppressed image.

[0044] The enhancement module is used to perform local contrast enhancement processing based on the artifact-suppressed image using a contrast-limited adaptive histogram equalization algorithm to obtain an enhanced image.

[0045] The calculation module is used to extract the three-dimensional gray-level histogram of the implanted target region from the enhanced image, and to calculate the gray-level skewness and gray-level kurtosis of the three-dimensional gray-level histogram.

[0046] The extraction module is used to extract features from the implantation target region in the enhanced image to obtain a texture feature set, wherein the texture feature set includes contrast, correlation, energy and homogeneity;

[0047] The analysis module is used to perform principal component analysis on the texture feature set to obtain principal component features that characterize the main variation direction of the trabecular bone structure. The gray level skewness, gray level kurtosis and the principal component features are concatenated according to preset weights to obtain a fused feature vector.

[0048] The verification module is used to analyze the fused feature vector using the density peak clustering algorithm to obtain the clustering results, map the clustering results to the anatomical partitions corresponding to the implantation target area, and perform Z-score verification on each sub-region corresponding to the clustering results in combination with the reference grayscale data to obtain the bone lesion area.

[0049] Thirdly, this application provides an electronic device, comprising:

[0050] Memory, used to store computer programs;

[0051] A processor is configured to execute the computer program to implement the steps of the target lesion identification method based on medical images as described in the first aspect above.

[0052] Fourthly, this application provides a computer-readable storage medium storing a computer program that, when executed by a processor, can implement the steps of the target lesion identification method based on medical images as described in the first aspect above.

[0053] The target lesion identification method based on medical imaging provided in this application has the following beneficial effects:

[0054] This application uses pre-implantation mandibular cone-beam computed tomography (CBCT) images and reference grayscale data to provide a baseline for subsequent analysis. A metal artifact reduction algorithm is then employed to effectively suppress implant artifacts, thereby improving image quality. Subsequently, local contrast enhancement based on the artifact-suppressed images improves the visual distinction between different bone regions. Next, the skewness and kurtosis of the three-dimensional grayscale histogram are extracted, and bone characteristics are quantified from a statistical distribution perspective. Simultaneously, a texture feature set including indicators such as contrast and correlation is extracted to characterize the microstructure of bone trabeculae. Principal component analysis of the texture features captures the main variation directions, and the above features are fused to construct a feature vector comprehensively describing bone quality. Finally, density peak clustering analysis, combined with statistical tests of anatomical partitions and reference data, enables objective and automated identification of lesion regions, improving the accuracy and stability of lesion identification.

[0055] Furthermore, after obtaining the initial region division through density peak clustering, this application sequentially performs statistical anomaly screening, boundary optimization combining anatomical rules and image gradients, calculates distribution deviation, and finally integrates multiple indicators to generate anomaly confidence. This process integrates spatial consistency, anatomical rationality, and statistical distribution characteristics, further improving the accuracy and robustness of lesion identification. Attached Figure Description

[0056] To more clearly illustrate the technical solutions of the embodiments of this application or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are only some embodiments of this application. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0057] Figure 1 A flowchart illustrating a target lesion identification method based on medical images, provided for an embodiment of this application;

[0058] Figure 2 A schematic diagram illustrating a specific implementation of a target lesion identification method based on medical images provided in this application embodiment;

[0059] Figure 3 A schematic diagram of the structure of a target lesion identification system based on medical images provided in this application embodiment;

[0060] Figure 4 This is a schematic diagram of the structure of an electronic device provided in an embodiment of this application. Detailed Implementation

[0061] In preoperative planning for dental implant surgery, accurate identification of lesion areas in the mandible is a key prerequisite for assessing bone quality and ensuring long-term implant stability. However, existing technologies have the following obvious shortcomings: on the one hand, visual interpretation, which relies on the doctor's experience, is highly subjective and has low repeatability; on the other hand, common automatic analysis methods usually rely on only a single image feature and have limited ability to distinguish the complex and variable microstructure of bone trabeculae and the interference of metal artifacts commonly found in images. This results in inaccurate and unreliable identification results.

[0062] To address this, this application proposes a target lesion identification method based on medical images. The technical concept of this method is as follows: First, targeted artifact suppression and local enhancement are applied to cone-beam CT images to improve image quality. Then, multiple features reflecting the overall statistical characteristics and microscopic texture of bone are extracted simultaneously and fused into a comprehensive feature description. Finally, a data-driven clustering algorithm combined with objective statistical tests based on anatomical location is used to achieve automated identification of lesion regions. Therefore, this scheme, by fusing multi-dimensional image information and introducing quantitative analysis, effectively improves the objectivity, accuracy, and adaptability to complex images in the identification process, thus solving the problems of subjectivity, bias, and weak anti-interference capabilities of existing methods.

[0063] To enable those skilled in the art to better understand the present application, the present application will be further described in detail below with reference to the accompanying drawings and specific embodiments. Obviously, the described embodiments are merely some embodiments of the present application, and not all embodiments. Based on the embodiments in this application, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of this application.

[0064] The core of this application is to provide a method for identifying target lesions based on medical images, and a flowchart of one specific implementation is shown below. Figure 1 As shown, the method includes:

[0065] S101. Obtain cone-beam CT images of the mandible before oral implant surgery and reference grayscale data corresponding to the implantation target area.

[0066] The implantation target area refers to the specific jawbone location where the implant is planned to be placed after clinical evaluation. This area is located within the mandible, and its precise spatial range can be determined by combining cone-beam CT images with oral examination results.

[0067] Reference grayscale data refers to a set of grayscale benchmarks for normal bone tissue used for quantitative comparison. Its typical source is the corresponding area of ​​the contralateral healthy mandible of the same patient, or a statistical model from a standard population database. This data is usually expressed in the form of grayscale mean and standard deviation of a specific area.

[0068] In step S101, firstly, a three-dimensional digital image of the patient's mandible is acquired using an oral cone-beam CT scanner. This image fully records the three-dimensional structure and density information of the mandible, including the implantation target area. Simultaneously, based on the implantation location determined by the established clinical plan, the three-dimensional spatial coordinates of the target area are automatically extracted from the image. Secondly, to establish a judgment benchmark, this application acquires reference grayscale data in parallel. This process is mainly accomplished by retrieving the grayscale information of the corresponding anatomical position of the contralateral healthy mandible of the same patient, or by querying a preset standard mandibular grayscale atlas database. It is necessary to ensure that the acquired reference data and the image of the target area are strictly aligned in anatomical position, thereby ensuring that subsequent comparative analysis is conducted under the same biological and imaging conditions.

[0069] S102. The cone-beam CT image is processed using a metal artifact reduction algorithm to obtain an artifact-suppressed image.

[0070] In one specific implementation, step S102 includes:

[0071] Step 1021: Read the original projection data of the cone-beam CT image to obtain the original sine wave, and perform metal region segmentation on the sine wave to obtain the metal region mask.

[0072] The original sine wave refers to the original two-dimensional projection data matrix collected by the detector during cone-beam CT scanning, arranged according to different rotation angles and detector unit positions. It records the attenuation information of X-rays after penetrating the object. The metal region mask is a binary image used to identify the projection area of ​​the metal object in the original sine wave. The pixel position where the metal projection is located is marked as the foreground, and the remaining area is marked as the background.

[0073] In step 1021, the unreconstructed original projection data is first read from the original data storage format of the cone-beam CT device. This data is arranged in order according to the scanning angle and detector channel and directly constitutes the original sine map. Subsequently, grayscale analysis is performed on the original sine map. Since the metal material has a strong attenuation ability for X-rays, it will form a local abnormally bright feature area on the original sine map. Then, by calculating the grayscale statistical characteristics of the entire sine map, a first segmentation threshold is dynamically set, and all pixels with grayscale values ​​higher than this threshold are determined to belong to the metal projection, thereby generating a binary matrix with the same size as the original sine map as a mask for the metal region.

[0074] For example, when processing mandibular scan data acquired by a certain type of cone-beam CT scanner, the original projection data is first read to obtain an original sine wave with a size of 600×300. Then, the maximum gray value Gmax of the sine wave is calculated, and then the first segmentation threshold is determined according to a preset rule. For example, the preset rule is: set the first segmentation threshold T to 80% of Gmax. For example, when Gmax is 3500, T=0.8×Gmax=0.8×3500=2800.

[0075] Next, each pixel of the sine wave is traversed, and pixels with a grayscale value greater than or equal to 2800 are marked as 1 in the mask matrix. Pixels with a value of 1 correspond to the areas in the original sine wave that are identified as metal projections. The remaining pixels are marked as 0, and pixels with a value of 0 represent non-metallic areas, i.e., projection areas of normal tissues such as bones and soft tissues. Finally, a 600×300 metal region mask is generated. The first segmentation threshold in the preset rules is an empirical value determined by analyzing a large amount of clinical data to effectively distinguish between metal and bone tissue. This application does not specifically limit the value of this threshold.

[0076] Step 1022: Based on the metal region mask, use the metal artifact reduction algorithm to extract multi-scale features from the original sine wave to obtain contextual features.

[0077] In step 1022, the explanation of the metal artifact reduction algorithm can be found in relevant technologies and will not be repeated here. The metal artifact reduction algorithm uses the metal region mask as key spatial guiding information to achieve targeted feature extraction and contextual modeling of the metal artifact region in the original sine wave image. Specifically:

[0078] The original sine wave and the metal region mask are fed into a neural network module containing multiple parallel convolutional paths. Each path uses a convolutional kernel of a different size. The metal region mask plays a role in guiding spatial attention in this process, enabling the network to focus on capturing subtle texture changes at the metal edges using small-sized convolutional kernels, while large-sized convolutional kernels are used to perceive the overall grayscale trend and structural patterns of a wider area.

[0079] The core task of this module is to analyze and learn the contextual features of normal tissue projection data around the metal projection area. These contextual features can include its texture regularity and spatial continuity. Then, the feature maps output by each path are spliced ​​and fused in the channel dimension to finally generate a feature tensor rich in multi-scale spatial contextual information.

[0080] For example, the original sine wave and the metal region mask are input into the multi-scale feature extraction module of the metal artifact reduction algorithm. This module consists of three parallel convolutional branches. The first branch uses a 3×3 convolutional kernel to extract local detail features; the second branch uses a 5×5 convolutional kernel to extract medium-range region features; and the third branch uses a 7×7 convolutional kernel to extract broader contextual features. After each branch performs convolution operations on the input data, a nonlinearity is introduced through the ReLU activation function. Subsequently, the feature maps output by the three branches are concatenated. Assuming each branch outputs 32 feature channels, the concatenated feature tensor generates a 96-channel feature tensor. The spatial size of this tensor is maintained by padding operations to preserve the same size as the original sine wave. Figure 1 This ensures that the spatial location information is fully preserved.

[0081] Step 1023: Input the feature sequence corresponding to the context feature and the projection statistical feature corresponding to the original sine curve into the gated recurrent unit for fusion to obtain the enhanced feature vector.

[0082] In step 1023, the projection statistical features corresponding to the original sine wave are first calculated. For example, the average value of the projection signal corresponding to each scanning angle can be calculated to form a projection statistical feature that describes the trend of data distribution along the angle. At the same time, the context features are reorganized into a feature sequence according to the order of the scanning angles. Then, the projection statistical features and the feature sequence are input into a gated recurrent unit network. The network processes the data under each scanning angle in sequence and uses its internal gating mechanism to dynamically fuse the context detail features and statistical trend features of the current angle, and remembers and transmits the dependencies between the sequences. The final output state is then further aggregated to form a comprehensive enhanced feature vector.

[0083] For example, for the original sine wave with a size of 600×300 in the example above, we first calculate the average gray value of each of the 300 scanning angle directions to obtain a statistical feature sequence containing 300 elements. At the same time, we convert the context feature tensor into a feature sequence of 300 elements according to the angle dimension. Then, we input these two sequences into a two-layer gated recurrent unit network. The network processes the input step by step and finally outputs an enhanced feature vector with a fixed dimension that integrates global statistical regularity and local context information.

[0084] Step 1024: Based on the enhanced feature vector, data reconstruction is performed using a conditional generative adversarial network to obtain the repair data corresponding to the metal region mask.

[0085] In step 1024, a pre-trained conditional generative adversarial network (GAN) model is used to complete the data repair. The generator of this model takes the obtained enhanced feature vector as the conditional input, and the network structure of the generator usually contains a series of deconvolution or upsampling layers to decode the input abstract feature vector and upsample it step by step, finally generating a two-dimensional data block that matches the size of the defect area identified by the metal area mask, i.e., the repair data, such as a 50×80 data matrix. It should be noted that the explanation of the conditional generative adversarial network can be found in related technologies, and its specific structural design, pre-training process and specific implementation process in practical applications are not described in this application.

[0086] Step 1025: Fuse the repaired data with the original data corresponding to the original sine curve to obtain an intermediate sine curve. Optimize the intermediate sine curve using an iterative gradient descent algorithm based on total variation minimization constraints to obtain the target sine curve.

[0087] Among them, the total variation minimization constraint is a commonly used image regularization constraint. Its core is to make the image smooth as a whole while allowing the preservation of necessary edge information.

[0088] In step 1025, the repair data is first filled into the corresponding area of ​​the original sine wave according to the precise coordinate position indicated by the metal area mask, replacing the original data marked as metal artifacts, thus forming a preliminary complete intermediate sine wave. Since direct replacement may cause discontinuities or abrupt transitions between the repaired area and the surrounding area, the intermediate sine wave needs to be optimized. That is, the gradient descent method is used to solve an optimization problem. The objective function of this problem aims to balance two objectives: first, the optimized sine wave should be as close as possible to the intermediate sine wave to retain information; second, the total variation of the optimized sine wave should be as small as possible to ensure overall smoothness. Then, through multiple iterations, a target sine wave is finally obtained that has a natural transition between the repaired area and the original area, and suppresses noise and inconsistencies.

[0089] For example, the generated repair data block is filled into the specified position of the intermediate sine wave; then, the smoothing constraint strength parameter q is set to 0.1, and the intermediate sine wave is used as the initial value. The gradient descent algorithm is used for iterative optimization. In each iteration, the gradient of the objective function is calculated and the pixel value of the sine wave is updated. After a set number of iterations, the algorithm converges and the optimized target sine wave is obtained. It should be noted that the explanation of the gradient descent algorithm can be found in related technologies, and its specific implementation process will not be described in this application.

[0090] Step 1026: Filter and back-project the target sine curve to obtain an artifact-suppressed image.

[0091] In step 1026, a commonly used filter, such as a ramp filter, can be selected to filter the target sine curve to perform frequency domain correction on the projection data before back projection. It should be noted that the explanation of the ramp filter can be referred to the relevant technology, and its specific implementation process will not be elaborated in this application.

[0092] The filtered projection data from each angle are then back-projected along their X-ray paths into a three-dimensional image space grid, and the back-projection results from all angles are accumulated to finally reconstruct a three-dimensional cone-beam CT body data image as an artifact-suppressed image. Since the metal artifacts in the input target sine wave have been effectively suppressed, the strip or star-shaped artifacts produced by metal are significantly reduced in the reconstructed three-dimensional image, and the skeletal structure is clearly presented.

[0093] This application suppresses metal artifacts in cone-beam CT images through the above process and obtains high-quality images with clear bone structure and significantly reduced artifacts, thus providing a reliable basis for subsequent accurate bone analysis and lesion identification.

[0094] S103. Based on the artifact-suppressed image, a contrast-limited adaptive histogram equalization algorithm is used to perform local contrast enhancement processing to obtain an enhanced image.

[0095] In one specific implementation, step S103 includes:

[0096] Step 1031: Perform regional characteristic analysis on the artifact suppression image to obtain a first characteristic region representing a high bone density region and a second characteristic region representing a low bone density region.

[0097] In step 1031, the global grayscale histogram of the artifact-suppressed image is first analyzed, and then a second segmentation threshold is automatically calculated based on the histogram to distinguish the two main grayscale distribution peaks. Then, the artifact-suppressed image is binarized and segmented using the threshold, and voxels with grayscale values ​​higher than the threshold are assigned to the first characteristic region, and voxels with grayscale values ​​lower than the threshold are assigned to the second characteristic region, thereby generating two independent region masks.

[0098] For example, for a mandibular artifact-suppressed image, its grayscale value ranges from 0 to 4095. Then, by analyzing its grayscale histogram, the second segmentation threshold R=2200 is calculated using the Otsu method. Then, all voxels with grayscale values ​​I>2200 in the image are marked as the first characteristic region, and voxels with grayscale values ​​I≤2200 are marked as the second characteristic region, where I represents the voxel grayscale value. It should be noted that the explanation of the Otsu method can be referred to the relevant technology, and its calculation process will not be described in detail in the embodiments of this application.

[0099] Step 1032: Generate regional statistical features to guide the global enhancement strategy based on the characteristic differences between the first characteristic region and the second characteristic region.

[0100] In step 1032, the mean and standard deviation of grayscale values ​​of all voxels within the first and second characteristic regions are calculated respectively. These two statistics describe the overall brightness and internal uniformity of regions with different bone density. The mean and standard deviation of the two regions are combined into a region statistical feature, which serves as the basis for setting parameters for subsequent enhancement algorithms.

[0101] For example, suppose that in the first characteristic region, the gray mean μ1 = 2800 and the standard deviation σ1 = 350 are calculated for all voxel points; in the second characteristic region, the gray mean μ2 = 1500 and the standard deviation σ2 = 250 are calculated. The mean and standard deviation of the two regions are combined to generate a vector [2800, 350, 1500, 250] representing the region's statistical features.

[0102] Step 1033: Based on the statistical characteristics of the region, apply a first constraint strategy, which mainly increases the clipping limit, to the first characteristic region using a contrast-limited adaptive histogram equalization algorithm, and apply a second constraint strategy, which mainly decreases the clipping limit, to the second characteristic region.

[0103] Among them, the contrast-limited adaptive histogram equalization algorithm is an image enhancement method that avoids excessive noise amplification by limiting the amplitude of local histogram equalization. Its cropping limit is a key parameter of the algorithm, and it sets the maximum number of pixels that a single gray level can accumulate in the local histogram.

[0104] In step 1033, the regional statistical features are transformed into specific adaptive histogram equalization algorithm parameter settings. Specifically, for the first characteristic region, a first clipping limit value is set based on the first standard deviation to suppress noise enhancement; for the second characteristic region, a second clipping limit value with an amplitude smaller than the first clipping limit value is set based on the standard deviation of the value being lower than the first standard deviation to improve the contrast enhancement effect. This parameter differentiation setting is achieved by directly adjusting the clipping limit parameters used by the adaptive histogram equalization algorithm when processing different regions.

[0105] For example, let the baseline cutting limit Cbase = 0.02, and the reference standard deviation be... =300, then calculate the clipping limit of the first characteristic region: And calculate the clipping limit of the second characteristic region: Where C1 and C2 are the first clipping limit value of the first characteristic region and the second clipping limit value of the second characteristic region, respectively. and These are the first standard deviation and the second standard deviation, respectively, for example. It is 350. It is 250.

[0106] Step 1034: Simultaneously, based on the preset anatomical reference images corresponding to the first characteristic region and the second characteristic region, the execution process of the first constraint strategy and the second constraint strategy is guided by a regularization term based on the gradient magnitude to obtain an intermediate image.

[0107] In step 1034, the artifact-suppressed image and the anatomical reference image are first registered and aligned in three-dimensional space. Then, the gradient magnitude maps of the corresponding regions of the two characteristic regions in the anatomical reference image are calculated respectively. During the enhancement process of parameterized adaptive histogram equalization for the two characteristic regions, gradient consistency constraints are applied simultaneously. That is, the difference between the gradient magnitude of the processed region image and the gradient magnitude of the corresponding region in the reference image is calculated as a regularization loss term. The enhancement intensity of the adaptive histogram equalization is then adjusted through the feedback of this loss term so that the enhancement process can improve the contrast while maintaining the original anatomical edges.

[0108] For example, the artifact-suppressed image is registered and aligned with the contralateral healthy mandibular CT image of the same patient, and the corresponding parts of the two characteristic regions in the anatomical reference image are extracted respectively. Then, the gradient magnitude maps Gref1 and Gref2 of the two characteristic regions in the anatomical reference image are calculated as reference benchmarks for the edge intensity of the normal anatomical structures in the corresponding regions. Then, after applying the first clipping limit value C1 to the first characteristic region for adaptive histogram equalization, the gradient magnitude map G1 of the processed region is calculated, and then the regularization loss term is calculated. Where N is the number of voxels in the corresponding region. If L1>0.001, then C1 is adjusted and reprocessed until L1≤0.001. Similarly, a similar process is performed on the second characteristic region. Through this regularization constraint based on gradient magnitude, while improving contrast, the edge intensity features that conform to normal anatomical structures are effectively maintained, and finally an intermediate image that satisfies the gradient consistency constraint is obtained.

[0109] Step 1035: Perform boundary smoothing processing on the intermediate image to obtain an enhanced image.

[0110] In step 1035, the three-dimensional boundary surface between the two regions is first identified. Then, a transition zone of fixed width is defined on both sides of the boundary surface, and the final gray value of the voxel in the transition zone is determined. The final gray value can be determined by linear interpolation of the gray values ​​after processing of the two regions. The interpolation weight depends on the distance of the voxel to the original boundary of the two regions. Then, spatial adaptive fusion is used to make the image transition smoothly at the boundary of the regions, thereby obtaining an enhanced image.

[0111] For example, suppose the width of the transition zone is 3 voxels, and for any voxel point P on the boundary surface, its minimum distance to the original boundary of the first characteristic region is d1, and its minimum distance to the original boundary of the second characteristic region is d2. This can be determined using the formula: Calculate the final grayscale value I1 and I2 are the values ​​of the intermediate image at the corresponding positions in the original first and second characteristic regions, respectively. Then, the calculation is completed by traversing all boundary voxels to enhance the image.

[0112] This application divides regions by analyzing differences in bone density, and adaptively adjusts the parameters of adaptive histogram equalization based on regional statistical characteristics to achieve targeted enhancement. Then, gradient constraints from anatomical reference images are introduced to maintain the authenticity of anatomical structures. Finally, boundary smoothing is applied to obtain enhanced images with prominent details and natural transitions to meet clinical diagnostic needs.

[0113] S104. Extract the three-dimensional grayscale histogram of the implantation target region from the enhanced image, and calculate the grayscale skewness and grayscale kurtosis of the three-dimensional grayscale histogram.

[0114] In one specific implementation, step S104 includes:

[0115] Step 1041: Locate the three-dimensional spatial range corresponding to the implantation target area in the enhanced image, and obtain the gray value of each voxel within the three-dimensional spatial range.

[0116] In step 1041, the range of the implantation target area is first precisely delineated in the enhanced image based on predefined three-dimensional spatial coordinates or through an image segmentation algorithm; then, the system traverses each voxel within the three-dimensional range and reads its corresponding digital grayscale value.

[0117] For example, based on the coordinate information output by the implantation planning software, a cubic region of 20mm×15mm×10mm is determined in the enhanced image as the implantation target region. For instance, this region contains a total of 60,000 voxel points. Subsequently, the system reads the grayscale value of each voxel point in sequence to obtain a set containing 60,000 grayscale values.

[0118] Step 1042: Based on all the gray values, a three-dimensional gray-level histogram of the implantation target region is statistically generated.

[0119] In step 1042, the number of times each possible gray value appears in the set is counted; then the gray value set is traversed and a count is established for each gray value that appears, and finally a frequency distribution map from the minimum gray value to the maximum gray value is formed as a three-dimensional gray histogram. This histogram intuitively reflects the central tendency and dispersion of bone gray values ​​in the target area.

[0120] For example, if the grayscale value range of the enhanced image is known to be 0 to 4095, then 60,000 grayscale values ​​are traversed and the number of occurrences of each grayscale level is counted. Suppose that after the count, the grayscale value 1500 appears 850 times, the grayscale value 1501 appears 820 times, and so on. Then, all grayscale levels and their corresponding frequencies are plotted to obtain the three-dimensional grayscale histogram of the implanted target area.

[0121] Step 1043: Based on the preset target grayscale range, the three-dimensional grayscale histogram is truncated and renormalized to obtain a corrected histogram. The preset target grayscale range is dynamically determined based on the mean and standard deviation of the reference grayscale data.

[0122] Here, truncation refers to forcibly setting the gray-level frequencies outside the interval in the histogram to zero; renormalization refers to rescaling the gray-level frequencies retained in the interval after truncation so that their sum is restored to 1, so as to ensure that the corrected histogram is still an effective probability distribution.

[0123] In step 1043, the target grayscale interval is first dynamically calculated based on the statistical characteristics of the reference grayscale data. Specifically, the lower limit of the interval... and upper limit Determined by the following formula: , ,in, It is the mean of the reference grayscale data. It is its standard deviation, and k is a preset positive coefficient, such as 3, used to control the range of interval coverage, which represents the main distribution range of gray values ​​in normal bone tissue; subsequently, the original three-dimensional gray-level histogram is truncated: all gray values ​​less than or greater than The frequency of the truncated gray level is set to 0; then, the gray level frequencies within the interval are renormalized to obtain the corrected histogram.

[0124] Step 1044: Calculate the second, third, and fourth central moments of the modified histogram.

[0125] Among them, the central moment is a statistic used in probability theory to describe the shape of the distribution; the second central moment, or variance, is used to describe the dispersion of the data distribution; the third central moment is related to skewness and is used to describe the asymmetry of the distribution; and the fourth central moment is related to kurtosis and is used to describe the sharpness of the distribution.

[0126] In step 1044, the mean of the distribution, i.e., the first central moment, is first calculated. Then, the second, third, and fourth central moments are calculated sequentially using the following formulas: , , , ,in, It is the value of the i-th gray level. It is its corresponding normalized frequency. , , These represent the second, third, and fourth order central moments, respectively, and the calculation process iterates through all gray levels within the target gray range.

[0127] Step 1045: Combine the third-order central moment and the second-order central moment to obtain the gray-level skewness; combine the fourth-order central moment and the second-order central moment to obtain the gray-level kurtosis.

[0128] Among them, gray skewness is a statistic that measures the asymmetry of gray distribution. Positive skewness indicates that there is a long tail on the right side of the distribution, while negative skewness indicates that there is a long tail on the left side. Gray kurtosis is a statistic that measures the steepness or flatness of the distribution curve. Kurtosis indicates that the distribution is more angular, while kurtosis indicates that the distribution is flatter.

[0129] In step 1045, the grayscale deviation is calculated using a standard formula. and grayscale kurtosis The calculation formula is as follows: , The calculated skewness and kurtosis are standardized dimensionless values, which facilitates subsequent unified quantitative analysis and comparison.

[0130] This application can sensitively reflect the symmetry and concentration changes in bone grayscale distribution, thus providing key features for subsequent quantitative assessment of bone homogeneity and health status.

[0131] S105. Extract features from the implantation target region in the enhanced image to obtain a texture feature set, wherein the texture feature set includes contrast, correlation, energy, and homogeneity.

[0132] In step S105, firstly, based on the three-dimensional range of the located implantation target area in the enhanced image, all voxel data within that area are extracted; then, for the three-dimensional image data, the statistical information of the texture is obtained by calculating its gray-level co-occurrence matrix; then, four core statistical features are extracted from the gray-level co-occurrence matrix: contrast, correlation, energy, and homogeneity. These four features together constitute a multi-dimensional texture description vector, i.e., a texture feature set.

[0133] This application can effectively characterize the microstructural properties of bone trabeculae, thus providing an important structural information basis for subsequent fusion of other features and lesion identification.

[0134] S106. Perform principal component analysis on the texture feature set to obtain principal component features that characterize the main variation direction of the trabecular structure. Concatenate the gray-level skewness, gray-level kurtosis and the principal component features according to preset weights to obtain a fused feature vector.

[0135] In one specific implementation, step S106 includes:

[0136] Step 1061: Decentralize the texture feature set to obtain a centralized feature set.

[0137] In step 1061, the average value of all samples in the texture feature set in each feature dimension is calculated; then, the corresponding feature mean is subtracted from each feature value of each sample in the original feature set to obtain a brand new centered feature set.

[0138] Step 1062: Based on the centered feature set, calculate the covariance matrix and obtain multiple eigenvectors and corresponding eigenvalues ​​of the covariance matrix.

[0139] In step 1062, after obtaining the centered feature set Xc, its covariance matrix is ​​calculated. For a dataset with n features, the size of the covariance matrix is ​​n×n. Subsequently, the covariance matrix is ​​decomposed into eigenvalues, which is a standard matrix operation that yields a set of paired eigenvalues ​​λi and eigenvectors vi.

[0140] Step 1063: Based on the variance contribution represented by the eigenvalues, select the top K eigenvectors with the highest variance contribution from the eigenvectors to obtain the principal component vectors.

[0141] Among them, the variance contribution refers to the proportion of a principal component eigenvalue to the sum of all eigenvalues. It reflects the amount of original data information carried by the principal component, and the cumulative variance contribution is the sum of the variance contributions of the first K principal components.

[0142] In step 1063, the variance contribution of each principal component is first calculated, and the contribution of the i-th principal component is... ,in, The sum of variances contributed by all principal components is calculated. Then, the eigenvalues ​​are sorted from largest to smallest and their corresponding variance contributions are accumulated sequentially. When the accumulated contribution first exceeds a preset threshold θ, such as θ=0.85, the first K eigenvectors v1, v2, ..., vK corresponding to this time are selected as principal component vectors to form the initial projection matrix VK.

[0143] Step 1064: Obtain preset anatomical reference data reflecting the trabeculae of the mandible, and based on the anatomical reference data, perform an element-wise multiplication weighted operation on the principal component vector to obtain an intermediate vector. The preset anatomical reference data is stored in the form of a standard trabecular orientation probability map.

[0144] Anatomical reference data refers to prior knowledge that provides the spatial arrangement of bone trabeculae. It is complementary to reference grayscale data. That is, reference grayscale data is used to describe the density statistical characteristics of normal bone, while anatomical reference data describes the spatial orientation probability of its microstructure.

[0145] In step 1064, a weight vector w corresponding to the anatomical position of the current implantation target area is first extracted from the standard trabecular bone direction probability map. The length of this vector is the same as the feature dimension, and the value of each element reflects the importance of the corresponding texture feature dimension in the normal anatomical structure. Furthermore, the elements in the weight vector w can be determined by normalizing the probability values ​​in the map. For example, the higher the probability of a certain direction, the greater the weight value of its corresponding feature dimension.

[0146] Subsequently, element-wise multiplication is performed on each principal component vector vi in ​​the initial projection matrix VK: This operation is equivalent to multiplying the initial projection matrix by a diagonal weight matrix W = diag(w) to obtain the weighted intermediate projection matrix. , as an intermediate vector.

[0147] Step 1065: Based on the intermediate vector, perform a projection transformation on the centered feature set to obtain a projected feature set, and then standardize the projected feature set to obtain principal component features.

[0148] In step 1065, the weighted intermediate projection matrix is ​​used. The linear projection transformation of the centered feature set Xc can be calculated using the formula: H = Xc × To obtain the dimension-reduced projected feature set H, since the scales of the new features after projection may be inconsistent, in order to eliminate the influence of dimensions, each feature of H can be standardized. Specifically, Z-score standardization is adopted, and finally the standardized principal component feature matrix Fpca with a mean of 0 and a standard deviation of 1 is obtained.

[0149] For example, the size of the centered feature set Xc is 100×4, and its weighted projection matrix is... The matrix H is 4×2 in size and is calculated by projection to be 100×2. Assuming that the mean of the first column of H is 0.05 and the standard deviation is 1.8, all elements in this column are transformed as follows: fnew=(fold−0.05) / 1.8, where fnew is the transformed element in the corresponding column and fold is the element before transformation in the corresponding column. A similar operation is then performed on the second column to finally obtain the standardized principal component feature Fpca.

[0150] This application extracts the main variation patterns of texture features through principal component analysis and innovatively introduces anatomical reference data to weight and guide the principal components, making the extracted features more consistent with the real biomechanical structure of bone trabeculae; ultimately, it constitutes a comprehensive descriptive vector that fully characterizes the distribution and microstructure of bone density, thereby improving the accuracy and robustness of subsequent lesion area identification.

[0151] S107. The density peak clustering algorithm is used to analyze the fused feature vector to obtain the clustering result. The clustering result is mapped to the anatomical partition corresponding to the implantation target area. The Z-score test is performed on each sub-region corresponding to the clustering result in combination with the reference grayscale data to obtain the bone lesion area.

[0152] In one specific implementation, such as Figure 2 As shown, step S107 includes:

[0153] Step 1071: Based on the fused feature vector, the density peak clustering algorithm is used to calculate the local density value and relative distance value of each data point.

[0154] The local density value refers to the number of other data points within a specified cutoff distance around each data point in the feature space, used to quantify the sample density of the neighborhood where the point is located; the relative distance value refers to the minimum distance from a data point to all data points with a higher local density value than it. For the point with the highest local density, its relative distance is defined as the maximum distance to all other points.

[0155] In step 1071, the fused feature vector set is used as input, where each fused feature vector corresponds to a data point in the implantation target region. For each data point, the algorithm calculates two core quantities: local density and local density. ,in, This represents the Euclidean distance between point i and point j. It is the preset cutoff distance. For indicator functions, in The function value is 1 when the condition is met, and 0 otherwise. Secondly, for point i with non-maximum local density, the function value is determined using the formula... Calculate the relative distance to this point That is, the minimum distance between point i and all points with higher density than it; otherwise, for point i with the maximum local density, it is determined by the formula... =max Calculate the relative distance to this point That is, the maximum distance from this point to all other points.

[0156] Step 1072: Combine the local density value and the relative distance value to construct a decision graph, identify multiple candidate cluster centers based on the decision graph, and obtain the clustering result based on the affiliation relationship between each data point and the candidate cluster centers.

[0157] The decision graph is a two-dimensional scatter plot, with the horizontal axis representing the local density values. The vertical axis represents the relative distance value. .

[0158] In step 1072, and This is plotted as a decision graph, in which potential cluster centers are distributed in... and In regions where the values ​​are all higher than the corresponding preset threshold, after determining the candidate cluster centers, each non-center point is assigned to the cluster center of the nearest point with a higher local density, thereby completing the clustering of the entire dataset.

[0159] Step 1073: Based on the clustering results, the implantation target region is divided into multiple clustered sub-regions, and each clustered sub-region is mapped to a corresponding anatomical partition. The Z-score is calculated by combining the mean gray value of the clustered sub-region with the mean gray value of the reference gray data.

[0160] In step 1073, a clustered subregion refers to the set of all voxel points that are classified into the same category in the three-dimensional image space; an anatomical partition is a spatial region predefined according to the standard mandibular anatomy; and the Z-score is used to quantify the degree of deviation of an observation from its reference population mean.

[0161] In step 1073, the cluster category labels are first mapped back to the 3D image space to form a 3D cluster sub-region. Then, by comparing the spatial coordinates, it is determined which anatomical region each cluster sub-region mainly falls into. Next, for each cluster sub-region, the average gray value of all voxels within it in the original enhanced image is calculated. Simultaneously, the mean value of the reference grayscale data corresponding to the same anatomical region is obtained. and standard deviation And according to the formula Calculate the Z-score value of this sub-region.

[0162] Step 1074: Based on the distribution of the local density value and the relative distance value corresponding to the clustered sub-regions, calculate the consistency metric within each clustered sub-region.

[0163] Among them, consistency measure is an index used to evaluate how closely members are distributed in the feature space within a clustered subregion.

[0164] In step 1074, for each clustered sub-region, the local density value set and relative distance value set of all data points belonging to that region are extracted; subsequently, the coefficient of variation of these two sets, i.e., the coefficient of variation of the local density values, are calculated respectively. ,in, and These are the standard deviation and mean of the local density value set, respectively. Similarly, the coefficient of variation of the relative distance values ​​is calculated. The coefficient of variation can be used to measure the degree of dispersion of members within a clustered subregion in terms of local density and relative distance.

[0165] The final consistency metric Con can be obtained using the following formula: Con ,in, and The preset weighting coefficients are, for example, respectively 0.5, 0.5, It is a very small positive number, such as 0.01.

[0166] Step 1075: Combine the preset anatomical constraint rules corresponding to the anatomical partition with the consistency metric, and introduce the structural tensor to optimize the boundaries of the corresponding clustered sub-regions to generate optimized sub-regions. The preset anatomical constraint rules are set based on the spatial position and morphology of the mandibular anatomical map, and the structural tensor is obtained based on image gradient calculation.

[0167] This application does not impose any restrictions on the pre-defined anatomical constraint rules, the content of the structural tensor, or its representation.

[0168] In step 1075, anatomical constraint rules are first applied to filter and correct each clustered sub-region. Then, a structure tensor is introduced for boundary optimization. For the sub-regions requiring optimization, the structure tensor field of their boundary voxels in the original enhanced image is calculated. The boundary optimization process is achieved by minimizing an energy function that encourages the optimized boundary to follow the strong image edges indicated by the structure tensor, while also considering the consistency metric Con. It should be noted that this embodiment does not limit the expression of this energy function and can be set accordingly based on actual conditions.

[0169] Step 1076: Calculate the distribution deviation between the grayscale distribution of the optimized sub-region and the normal grayscale model in the reference grayscale data.

[0170] In step 1076, firstly, the gray values ​​of all voxels within each optimized sub-region are statistically analyzed to form the empirical gray-level distribution of that region, denoted as distribution A(g), representing the probability of a voxel with a gray value of g appearing in the optimized sub-region. Simultaneously, a parameterized normal gray-level model B(g) is provided for each anatomical partition based on the gray-level data. This model is typically assumed to follow a Gaussian distribution, with the probability density function: B(g) ∼ N(μB, σ²). ), where μB and σ These are the mean and standard deviation of the grayscale values ​​for the anatomical region in the reference grayscale data, respectively. Subsequently, the Kullback-Leibler divergence is used as the method for calculating the distribution deviation Div. This calculation measures the difference between distribution A(g) and the normal grayscale model B(g), and the larger the Div value, the greater the difference between the grayscale distribution of the optimized sub-region and the normal model.

[0171] Step 1077: The Z-score, the consistency metric, and the corresponding distribution deviation are weighted and summed to obtain the anomaly confidence of each optimized sub-region.

[0172] Step 1077 may specifically include the following steps:

[0173] Step a1: Input the Z-score, the consistency metric, and the distribution deviation into a pre-trained three-layer fully connected neural network, and predict the fusion coefficients of each item through the three-layer fully connected neural network.

[0174] The structural design, parameter design, and training process of the three-layer fully connected neural network are not the focus of this application, and therefore will not be elaborated upon.

[0175] In step a1, a pre-trained three-layer fully connected neural network is used to learn the optimal fusion weights for the three input metrics. The input layer of the network receives three normalized scalars: the Z-score z, the consistency measure Con, and the distribution deviation Div. The output layer of the network uses the Softmax activation function to ensure that the three output fusion coefficients a, b, and c are non-negative and their sum is 1.

[0176] Step a2: Based on the fusion coefficient, the Z-score value, the consistency metric, and the distribution deviation are weighted and summed to obtain a comprehensive anomaly index.

[0177] In step a2, after obtaining the fusion coefficients, a weighted sum is calculated to obtain the comprehensive anomaly index, Score. .

[0178] Step a3: Convert the composite anomaly index into an intermediate confidence level using the Sigmoid function.

[0179] In step a3, to transform the comprehensive anomaly index Score into an intermediate value that conforms to a probabilistic interpretation, this step can apply the Sigmoid function for nonlinear transformation, the specific formula of which is: ,in, is the median confidence level, m is the slope parameter, and n is the bias parameter.

[0180] Step a4: Based on the prior Gaussian distribution of normal gray values ​​in the reference gray data, the posterior probability of the intermediate confidence is calculated using Bayes' formula to obtain the abnormal confidence of each optimized sub-region.

[0181] In step a4, Bayesian inference is introduced to determine the intermediate confidence level. This is treated as a likelihood estimate and corrected using prior clinical knowledge to obtain the posterior confidence level for abnormality. Let be the anomaly confidence level of the corresponding optimized sub-region, and let the prior probability of finding an anomaly region in this anatomical partition be . The prior probability of the normal region is =1- Its Bayesian calculation formula is: .

[0182] Step 1078: Select candidate lesion regions whose abnormal confidence exceeds a preset abnormal threshold, and spatially aggregate the candidate lesion regions to obtain bone lesion regions.

[0183] In step 1078, a global anomaly confidence threshold is first set. Traverse all optimized sub-regions and calculate their posterior anomaly confidence. Compare with a threshold, ≥ The regions are marked as candidate lesion regions; then, spatial aggregation is performed to form meaningful clinical lesions, and the three-dimensional spatial centroids of all candidate regions are calculated. If the Euclidean distance between the centroids of any two regions is less than the preset aggregation distance Δ, the two regions are merged into a new connected region, and this process is iterated until there are no more regions that can be merged; finally, all independent or merged connected three-dimensional regions are output as the final bone lesion regions.

[0184] For example, suppose there are five optimization sub-regions, and their posterior anomaly confidence levels are... The values ​​were 0.109, 0.803, 0.721, 0.450, and 0.915, respectively, and then a threshold was set. =0.7, and then compare, the three regions with posterior abnormal confidence scores of 0.803, 0.721, and 0.915 were selected as candidate lesion regions; then the aggregation distance Δ = 5.0 mm was set, and it was found that the centroid distance between the two regions with confidence scores of 0.803 and 0.721 was 3.2 mm, which was less than Δ, so they were merged; the region with confidence score of 0.915 was more than Δ away from the other regions and remained independent; finally, the system output two bone lesion regions: one was formed by merging the regions corresponding to posterior abnormal confidence scores of 0.803 and 0.721, and the other was an independent region with confidence score of 0.915.

[0185] This application quantifies regional anomalies from a statistical distribution perspective by calculating distribution deviation, and designs a decision fusion framework based on neural networks and Bayesian inference. It then fuses multi-source information into a robust anomaly confidence level. Finally, it outputs continuous lesion regions that conform to clinical cognition through threshold screening and spatial aggregation, thereby improving the accuracy and interpretability of lesion identification.

[0186] Figure 3 This is a schematic diagram illustrating a specific implementation of a target lesion identification system based on medical images, as provided in this application. Figure 3 The system may include:

[0187] The acquisition module 31 is used to acquire cone-beam CT images of the mandible before oral implant surgery and reference grayscale data corresponding to the implantation target area.

[0188] The processing module 32 is used to process the cone-beam CT image using a metal artifact reduction algorithm to obtain an artifact-suppressed image.

[0189] Enhancement module 33 is used to perform local contrast enhancement processing based on the artifact-suppressed image using a contrast-limited adaptive histogram equalization algorithm to obtain an enhanced image.

[0190] The calculation module 34 is used to extract the three-dimensional grayscale histogram of the implantation target region from the enhanced image, and to calculate the grayscale skewness and grayscale kurtosis of the three-dimensional grayscale histogram.

[0191] The extraction module 35 is used to extract features from the implantation target region in the enhanced image to obtain a texture feature set, wherein the texture feature set includes contrast, correlation, energy and homogeneity.

[0192] Analysis module 36 is used to perform principal component analysis on the texture feature set to obtain principal component features that characterize the main variation direction of the trabecular structure. The gray-level skewness, gray-level kurtosis and the principal component features are concatenated according to preset weights to obtain a fused feature vector.

[0193] The verification module 37 is used to analyze the fused feature vector using the density peak clustering algorithm to obtain the clustering result, map the clustering result to the anatomical partition corresponding to the implantation target area, and perform Z-score test on each sub-region corresponding to the clustering result in combination with the reference grayscale data to obtain the bone lesion area.

[0194] The target lesion identification system based on medical images in this application is used to implement the aforementioned target lesion identification method based on medical images. Therefore, the specific implementation of the target lesion identification system based on medical images can be found in the embodiment section of the target lesion identification method based on medical images above. The specific implementation can be referred to the description of the corresponding embodiments, which will not be repeated here.

[0195] like Figure 4 As shown, this application also provides an electronic device, including: a memory 41 for storing a computer program; and a processor 42 for executing the computer program to implement the steps of any of the above-described target lesion identification methods based on medical images.

[0196] This application also provides a computer-readable storage medium storing a computer program, which, when executed by a processor, implements the steps of any of the above-described target lesion identification methods based on medical images.

[0197] In one exemplary embodiment, the aforementioned computer-readable storage medium may include, but is not limited to, various media capable of storing computer programs, such as USB flash drives, read-only memory, random access memory, portable hard drives, magnetic disks, or optical disks.

[0198] Embodiments of the present invention also provide a computer program product, which includes a computer program that, when executed by a processor, implements the steps in any of the above embodiments of the target lesion identification method based on medical images.

[0199] Those skilled in the art will further recognize that the units and algorithm steps of the various examples described in conjunction with the embodiments disclosed herein can be implemented in electronic hardware, computer software, or a combination of both. To clearly illustrate the interchangeability of hardware and software, the components and steps of the various examples have been generally described in terms of functionality in the foregoing description. Whether these functions are implemented in hardware or software depends on the specific application and design constraints of the technical solution. Those skilled in the art can use different methods to implement the described functions for each specific application, but such implementations should not be considered beyond the scope of this invention.

[0200] The above provides a detailed description of the target lesion identification method and system based on medical images provided in this application. Specific examples have been used to illustrate the principles and implementation methods of this application. The descriptions of the embodiments above are only for the purpose of helping to understand the method and its core ideas. It should be noted that those skilled in the art can make several improvements and modifications to this application without departing from the principles of this application, and these improvements and modifications also fall within the protection scope of this application.

Claims

1. A method for identifying target lesions based on medical images, characterized in that, include: Obtain cone-beam CT images of the mandible before oral implant surgery and reference grayscale data corresponding to the implantation target area; The cone-beam CT images are processed using a metal artifact reduction algorithm to obtain artifact-suppressed images; Based on the artifact-suppressed image, a contrast-limited adaptive histogram equalization algorithm is used to perform local contrast enhancement processing to obtain an enhanced image. Extract the three-dimensional grayscale histogram of the implantation target region from the enhanced image, and calculate the grayscale skewness and grayscale kurtosis of the three-dimensional grayscale histogram; Feature extraction is performed on the implantation target region in the enhanced image to obtain a texture feature set, wherein the texture feature set includes contrast, correlation, energy and homogeneity; Principal component analysis is performed on the texture feature set to obtain principal component features that characterize the main variation direction of the trabecular bone structure. The gray-level skewness, gray-level kurtosis and the principal component features are concatenated according to preset weights to obtain a fused feature vector. The density peak clustering algorithm is used to analyze the fused feature vector to obtain the clustering results. The clustering results are mapped to the anatomical partitions corresponding to the implantation target area. The Z-score test is performed on each sub-region corresponding to the clustering results in combination with the reference grayscale data to obtain the bone lesion area. The density peak clustering algorithm is used to analyze the fused feature vector to obtain clustering results. These clustering results are then mapped to the anatomical regions corresponding to the implantation target area. A Z-score test is performed on each sub-region corresponding to the clustering results using the reference grayscale data to obtain the bone lesion region, including: Based on the fused feature vector, the density peak clustering algorithm is used to calculate the local density value and relative distance value of each data point; By combining the local density value and the relative distance value, a decision map is constructed. Based on the decision map, multiple candidate cluster centers are identified, and the clustering result is obtained according to the affiliation relationship between each data point and the candidate cluster centers. Based on the clustering results, the implantation target region is divided into multiple clustered sub-regions, and each clustered sub-region is mapped to a corresponding anatomical partition. The Z-score is calculated by combining the mean gray value of the clustered sub-region with the mean gray value of the reference gray data. Based on the distribution of the local density values ​​and the relative distance values ​​corresponding to the clustered sub-regions, a consistency metric is calculated within each clustered sub-region. Combining the preset anatomical constraint rules corresponding to the anatomical partitions with the consistency metric, and introducing the structural tensor to optimize the boundaries of the corresponding clustered sub-regions, an optimized sub-region is generated. The preset anatomical constraint rules are based on the spatial location and morphology of the mandibular anatomical diagram, and the structural tensor is obtained based on image gradient calculation. Calculate the distribution deviation between the grayscale distribution of the optimized sub-region and the normal grayscale model in the reference grayscale data; The Z-score, the consistency metric, and the corresponding distribution deviation are weighted and summed to obtain the anomaly confidence of each optimized sub-region. Candidate lesion regions with an abnormal confidence level exceeding a preset abnormal threshold are selected, and spatial aggregation is performed on the candidate lesion regions to obtain bone lesion regions.

2. The method for identifying target lesions based on medical images according to claim 1, characterized in that, The step of weighted summing of the Z-score, the consistency metric, and the corresponding distribution deviation to obtain the anomaly confidence of each optimized sub-region includes: The Z-score, the consistency metric, and the distribution deviation are input into a pre-trained three-layer fully connected neural network, which then predicts the fusion coefficients of each item. Based on the fusion coefficient, the Z-score, the consistency metric, and the distribution deviation are weighted and summed to obtain a comprehensive anomaly index. The composite anomaly index is converted into an intermediate confidence level using the Sigmoid function; Based on the prior Gaussian distribution of normal gray values ​​in the reference gray data, the Bayesian formula is used to calculate the posterior probability of the intermediate confidence level to obtain the abnormal confidence level of each optimized sub-region.

3. The method for identifying target lesions based on medical images according to claim 1, characterized in that, The cone-beam CT image is processed using a metal artifact reduction algorithm to obtain an artifact-suppressed image, including: The original projection data of the cone-beam CT image is read to obtain the original sine wave. The sine wave is then segmented into metal regions to obtain a metal region mask. Based on the metal region mask, a metal artifact reduction algorithm is used to extract multi-scale features from the original sine wave to obtain contextual features; The feature sequence corresponding to the context features and the projection statistical features corresponding to the original sine curve are input into a gated recurrent unit for fusion to obtain an enhanced feature vector. Based on the enhanced feature vector, data reconstruction is performed using a conditional generative adversarial network to obtain the repair data corresponding to the metal region mask. The repaired data is fused with the original data corresponding to the original sine curve to obtain an intermediate sine curve. The intermediate sine curve is then optimized using a gradient descent iterative algorithm based on total variation minimization constraints to obtain the target sine curve. The target sine wave is filtered and back-projected to reconstruct the image, resulting in an artifact-suppressed image.

4. The method for identifying target lesions based on medical images according to claim 1, characterized in that, The image based on the artifact-suppressed image is then subjected to local contrast enhancement processing using a contrast-limited adaptive histogram equalization algorithm to obtain an enhanced image, including: Regional characteristic analysis was performed on the artifact-suppressed image to obtain a first characteristic region representing a high bone density region and a second characteristic region representing a low bone density region. Based on the characteristic differences between the first characteristic region and the second characteristic region, regional statistical features are generated to guide the global enhancement strategy. Based on the statistical characteristics of the region, a first constraint strategy, which mainly increases the clipping limit, is applied to the first characteristic region through a contrast-limited adaptive histogram equalization algorithm, and a second constraint strategy, which mainly decreases the clipping limit, is applied to the second characteristic region. Meanwhile, based on the preset anatomical reference images corresponding to the first and second characteristic regions respectively, the execution process of the first and second constraint strategies is guided by a regularization term based on the gradient magnitude to obtain intermediate images. The intermediate image is then subjected to boundary smoothing processing to obtain an enhanced image.

5. The method for identifying target lesions based on medical images according to claim 1, characterized in that, The principal component analysis of the texture feature set yields principal component features characterizing the main variation directions of the trabecular bone structure, including: The texture feature set is decentralized to obtain a centralized feature set; Based on the centered feature set, the covariance matrix is ​​calculated, and multiple eigenvectors and corresponding eigenvalues ​​of the covariance matrix are obtained. Based on the variance contribution represented by the eigenvalues, the top K eigenvectors with the highest variance contribution are selected from the eigenvectors to obtain the principal component vectors; Obtain preset anatomical reference data reflecting the trabeculae of the mandible, and based on the anatomical reference data, perform an element-wise multiplication weighted operation on the principal component vector to obtain an intermediate vector. The preset anatomical reference data is stored in the form of a standard trabecular orientation probability map. Based on the intermediate vector, the centered feature set is projected to obtain the projected feature set, and the projected feature set is standardized to obtain the principal component features.

6. The method for identifying target lesions based on medical images according to claim 1, characterized in that, The step of extracting the three-dimensional grayscale histogram of the implantation target region from the enhanced image and calculating the grayscale skewness and grayscale kurtosis of the three-dimensional grayscale histogram includes: Locate the three-dimensional spatial range corresponding to the implantation target region in the enhanced image, and obtain the gray value of each voxel point within the three-dimensional spatial range; Based on all the gray values, a three-dimensional gray-level histogram of the implantation target region is statistically generated. Based on a preset target grayscale range, the three-dimensional grayscale histogram is truncated and renormalized to obtain a corrected histogram. The preset target grayscale range is dynamically determined based on the mean and standard deviation of the reference grayscale data. Calculate the second, third, and fourth central moments of the modified histogram; By combining the third-order central moment and the second-order central moment, the gray-level skewness is obtained; by combining the fourth-order central moment and the second-order central moment, the gray-level kurtosis is obtained.

7. A target lesion identification system based on medical images, characterized in that, A method for identifying target lesions based on medical images as described in any one of claims 1 to 6 includes: The acquisition module is used to acquire cone-beam CT images of the mandible before oral implant surgery and reference grayscale data corresponding to the implantation target area; The processing module is used to process the cone-beam CT image using a metal artifact reduction algorithm to obtain an artifact-suppressed image. The enhancement module is used to perform local contrast enhancement processing based on the artifact-suppressed image using a contrast-limited adaptive histogram equalization algorithm to obtain an enhanced image. The calculation module is used to extract the three-dimensional gray-level histogram of the implanted target region from the enhanced image, and to calculate the gray-level skewness and gray-level kurtosis of the three-dimensional gray-level histogram. The extraction module is used to extract features from the implantation target region in the enhanced image to obtain a texture feature set, wherein the texture feature set includes contrast, correlation, energy and homogeneity; The analysis module is used to perform principal component analysis on the texture feature set to obtain principal component features that characterize the main variation direction of the trabecular bone structure. The gray level skewness, gray level kurtosis and the principal component features are concatenated according to preset weights to obtain a fused feature vector. The verification module is used to analyze the fused feature vector using the density peak clustering algorithm to obtain the clustering results, map the clustering results to the anatomical partitions corresponding to the implantation target area, and perform Z-score verification on each sub-region corresponding to the clustering results in combination with the reference grayscale data to obtain the bone lesion area.

8. An electronic device, characterized in that, include: Memory, used to store computer programs; A processor, configured to execute the computer program to implement the steps of the target lesion identification method based on medical images as described in any one of claims 1 to 6.

9. A computer-readable storage medium, characterized in that, The computer-readable storage medium stores a computer program that, when executed by a processor, enables the implementation of the target lesion identification method based on medical images as described in any one of claims 1 to 6.

Images