A medical image processing method and system fusing imageomics and machine learning

Abstract

The application relates to the technical fields of medical image processing and artificial intelligence, and discloses a medical image processing method and system fusing imageomics and machine learning. The method comprises the following steps: acquiring pre-stored three-dimensional medical images, registration masks and structured attribute data, determining and solidifying resampling, gray limit and discretization parameters according to training samples; extracting statistical, shape and texture features through an original image, a spatial filter and a multi-scale transformation channel, sequentially applying stability, sparsity and correlation constraints to form a target feature subspace; fusing target image features or weighted scores thereof and structured attribute feature codes, inputting a classification model to obtain a classification reference value; calculating feature contribution degrees and interaction intensities, generating an interpretable data package containing the classification reference value, the contribution degree sequence, interaction data, a model version and pre-processing parameters, and outputting the same. The application is suitable for computer-aided analysis of pre-acquired three-dimensional medical images.

Description

Technical Field

[0001] This invention relates to the fields of medical image processing and artificial intelligence, specifically to a method, system, electronic device, and computer-readable storage medium for standardizing, extracting quantitative features through multiple pathways, filtering cascaded features, fusing structured data, and outputting interpretable classifications from pre-acquired three-dimensional medical images. Background Technology

[0002] Abdominal computed tomography (CT) images typically contain multidimensional information such as pixel intensity, morphology, and spatial texture. Radiomics techniques can convert three-dimensional voxel data of a region of interest into quantitative features and establish a mapping relationship between features and predetermined image categories through machine learning models. Taking a gastric lesion image that does not show clear signs of peritoneal abnormalities on conventional CT as an example, the primary lesion may still contain microscopic textures and spatial heterogeneity information that are not easily visible to the naked eye, and therefore can be used as an object for computer image analysis.

[0003] Existing radiomics processing schemes still have several technical problems. First, the slice thickness, voxel spacing, tube current, and reconstruction settings may vary between scans. Directly extracting texture features can easily lead to non-targeted fluctuations caused by acquisition parameters. Second, a single region of interest can generate hundreds or even thousands of high-dimensional features, while the number of training samples is usually limited. If the normalization parameters are fitted to all data, only a single screening algorithm is used, or the correlation between features is not controlled, data leakage, redundant computation, and model instability may be introduced. Third, while single linear models are easy to interpret, their ability to express nonlinear thresholds and feature interactions is limited; while nonlinear ensemble models can fit complex relationships, their outputs are not easily traceable. Fourth, some schemes only output high and low class labels, without simultaneously providing feature contribution, inter-model differences, and interaction strength, which is not conducive to verifying the algorithm's behavior.

[0004] Furthermore, existing technologies already employ enhanced CT images, radiomics features, and regression scores for abdominal image classification. Therefore, it remains necessary to provide a processing scheme with clearly defined technical constraints regarding image standardization, feature selection order, and interpretable output structure, enabling functional cooperation between algorithmic features and medical image processing features, and allowing for implementation by those skilled in the art.

[0005] The objects processed by this invention are digital images and structured data that have been acquired and stored through an independent process. The method itself does not involve scanning, sampling, surgery, medication or other diagnostic and treatment operations on living human or animal bodies. The output is image classification reference data for subsequent analysis, and does not directly provide disease diagnosis conclusions or treatment plans. Summary of the Invention

[0006] The technical problem that this invention aims to solve is: how to obtain a relatively stable and low-redundancy fusion feature expression when there are differences in the acquisition parameters of three-dimensional medical images and the candidate feature dimension is high, and to provide traceable feature contribution and interaction information while generating classification reference values.

[0007] To address the aforementioned technical problems, this invention provides a medical image processing method that integrates radiomics and machine learning. This method is executed by a computer on pre-acquired and stored three-dimensional medical image data, and includes data reading, image standardization, parallel feature extraction, cascaded feature filtering, structured feature fusion, classification mapping, and data interpretation generation.

[0008] During the image standardization stage, the 3D medical images are resampled to a preset isotropic voxel spacing, grayscale values ​​are limited, and discretization is performed using a fixed bin width or a fixed number of bins. The isotropic voxel spacing, grayscale upper and lower limits, and discretization parameters can be predetermined based on the imaging modality, voxel spacing distribution, and grayscale distribution of the training samples, and are saved as a configuration file associated with the model version. The application stage directly calls the configuration file to maintain consistency between the preprocessing conditions of the training and application stages.

[0009] In the parallel feature extraction stage, the computational domain is defined by a region of interest mask registered with the medical image space, and at least the original image branch, spatial filtering branch, and multi-scale transformation branch are constructed. Each branch outputs first-order statistical, shape, and texture features, respectively, and forms an image omics feature vector with a fixed field order according to a preset feature dictionary.

[0010] In the cascaded feature selection stage, the feature vector is first transformed using normalization parameters determined solely by the training samples. Then, stability constraints, sparsity constraints, and correlation constraints are applied sequentially. Stability constraints can be based on intra-group correlation coefficients obtained from repeated segmentation or multi-observer segmentation; sparsity constraints can be based on near-zero variance detection and regression with regularization terms; and correlation constraints can be based on Spearman or Pearson correlation coefficients. The stability threshold, regularization parameter, and correlation threshold can all be determined within the training samples, and the final dimensionality can be limited based on the ratio of positive events to predictor variables to ensure that the verification samples and application data do not participate in the refitting of these parameters.

[0011] In the fusion and classification stage, the target feature subspace may contain at least two target image features determined from the candidate feature set through cascaded screening. These target image features can be directly encoded as vectors, or a radiomics score can be formed according to R = β0 + β1 × f1 + ... + βm × fm, where m is a positive integer not less than 2, f1 to fm are normalized target image features, and β0 to βm are fitted from training samples and fixed with each model version. The target image feature vector or radiomics score is encoded and concatenated with legally and compliantly obtained structured attribute features, and then fed into at least one trained classification model. The classification model may include linear classification models and nonlinear ensemble classification models, and the main output branch can be determined by the discrimination, calibration, and precision-recall metrics obtained through cross-validation. The output of the classification model is represented as the degree of matching for a predetermined image category or a classification reference value.

[0012] During the interpretation phase, the contribution of each feature component is calculated using training samples or a representative subset of training samples as background data. For nonlinear models, the feature interaction strength is further calculated based on a bivariate partial dependency function. The system organizes classification reference values, feature contributions, contribution ranking, interaction matrix, model version, and preprocessing parameters into an interpretable data package, which can be written to local storage, a database, an image archiving communication system, or a display terminal.

[0013] The present invention also provides a medical image processing system, electronic device and non-transitory computer-readable storage medium corresponding to the above method.

[0014] Compared to schemes that only perform general imagemics scoring, the technical solution of this invention can produce the following beneficial effects: First, isotropic resampling, grayscale limiting, and fixed discretization parameters help reduce feature fluctuations caused by different layer thicknesses and intensity ranges; Second, the sequential combination of stability, sparsity, and correlation constraints can reduce redundant features and subsequent model computations while controlling data leakage; Third, the original image, spatial filtering, and multi-scale transformation branches can extract complementary information from different spatial frequency ranges; Fourth, the fusion of target image features and structured attribute features can provide feature representations different from those from a single data source; Fifth, linear and nonlinear models are interpreted in the same input space, which can identify differences in model outputs and non-additive interactions; Sixth, interpretable data packets retain model versions, standardized parameters, and contributing data, which helps with result verification, auditing, and subsequent model maintenance. The above effects are affected by the scale of training data, image quality, and specific model parameters, and should not be construed as an absolute guarantee of performance in any application scenario. Attached Figure Description

[0015] Figure 1 This is a schematic diagram of the overall functional modules of the medical image processing system of the present invention;

[0016] Figure 2 This is a schematic diagram of the overall flow of the medical image processing method of the present invention;

[0017] Figure 3 This is a schematic diagram of the image standardization and multi-path feature extraction process of the present invention;

[0018] Figure 4 This is a schematic diagram of the cascaded feature selection and model training process of the present invention;

[0019] Figure 5 This is a schematic diagram of the linear and nonlinear dual-path classification and interpretation output process of the present invention;

[0020] Figure 6 This is a schematic diagram of the hardware structure of the electronic device of the present invention.

[0021] Figure labeling: 10—Data receiving module; 20—Image preprocessing module; 30—Region of interest management module; 40—Feature extraction module; 41—Original image branch; 42—Spatial filtering branch; 43—Multi-scale transformation branch; 50—Feature filtering module; 51—Stability filtering unit; 52—Sparseness filtering unit; 53—Correlation filtering unit; 60—Fusion classification module; 61—Linear classification branch; 62—Nonlinear ensemble classification branch; 70—Interpretation module; 71—Contribution calculation unit; 72—Interaction strength calculation unit; 73—Data packet generation unit; 80—Output module; 90—Electronic device; 91—Processor; 92—Memory; 93—Medical image data interface; 94—Display interface; 95—Network interface; 96—Bus. Detailed Implementation

[0022] The present invention will be further described below with reference to the accompanying drawings and embodiments. The embodiments are used to illustrate the technical concept of the present invention and are not intended to unduly limit the scope of protection. Unless otherwise stated, the thresholds, model types, and software versions in the embodiments can be adaptively adjusted while meeting the same data processing purpose. Example 1: System Overall Structure

[0023] like Figure 1 As shown, the medical image processing system includes a data receiving module 10, an image preprocessing module 20, a region of interest management module 30, a feature extraction module 40, a feature filtering module 50, a fusion classification module 60, an interpretation module 70, and an output module 80. These modules can be deployed on the same computer device or distributed according to a client-server architecture.

[0024] The data receiving module 10 reads pre-acquired 3D medical image data, region-of-interest masks, and structured attribute data via a DICOM interface, file interface, or database interface. During reading, the data receiving module 10 verifies image sequence identifiers, voxel spacing, layer number, grayscale units, mask size, and spatial transformation matrix, and writes the verification results to the task log. For data that does not meet the integrity requirements, the system can return a pending review status without generating classification reference values.

[0025] The image preprocessing module 20 performs resampling, grayscale limiting, and discretization according to the configuration file bound to the model version. The region of interest management module 30 resamples the mask to the normalized image coordinate system and checks whether the mask is empty, whether it exceeds the image boundary, and whether the number of voxels in the connected region reaches a preset lower limit.

[0026] The feature extraction module 40 includes an original image branch 41, a spatial filtering branch 42, and a multi-scale transformation branch 43. These three branches can run in parallel, and the output is written to a unified feature dictionary according to feature names and branch identifiers. The feature filtering module 50 sequentially calls the stability filtering unit 51, the sparsity filtering unit 52, and the correlation filtering unit 53. The fusion classification module 60 concatenates the target image features with structured attribute features and then performs numerical mapping. The interpretation module 70 includes a contribution calculation unit 71, an interaction strength calculation unit 72, and a data packet generation unit 73. The output module 80 is used to persist or visualize the data packets.

[0027] like Figure 2 As shown, the medical image processing method of the present invention is executed by a computer device on pre-acquired and stored data, including: reading three-dimensional medical image data, a region of interest mask spatially registered with the three-dimensional medical image data, and structured attribute data; performing voxel resampling, grayscale limiting, and grayscale discretization on the three-dimensional medical image data; extracting statistical features, shape features, and texture features through original image branches, spatial filtering branches, and multi-scale transformation branches; transforming the radiomics features using normalization parameters determined only by training samples, and sequentially applying stability constraints, sparsity constraints, and correlation constraints to form a target feature subspace; encoding and concatenating the target image features and structured attribute features to form a fused feature vector; inputting the fused feature vector into a trained classification model to obtain classification reference values; calculating the contribution and interaction strength of each feature component, generating an interpretable data package containing classification reference values, contribution ranking, interaction data, and model version information, and outputting the interpretable data package to a storage device or display device. The above processing steps are further described in subsequent embodiments. Example 2: Standardization of 3D Images

[0028] This embodiment uses pre-acquired and stored abdominal contrast-enhanced CT delayed-phase images. Exemplary acquisition parameters are: tube voltage 120 kV, automatic tube current modulation range 50–413 mAs, rotation time 0.5 s, matrix 512×512, original slice thickness 0.625–5 mm, and delayed-phase images acquired approximately 120–150 s after contrast agent injection. The above acquisition actions are not part of the computer processing method of this invention. Images obtained from other devices or with other parameters can also be used as input after configuration verification.

[0029] The image preprocessing module 20 reads the DICOM grayscale mapping parameters and converts the pixel values ​​into CT values. In this embodiment, trilinear interpolation is used to resample the image to 1 mm × 1 mm × 1 mm isotropic voxels. Nearest neighbor interpolation is used for the discrete label mask to avoid generating non-integer labels. The image voxel values ​​are limited to the range of -1000 HU to 2000 HU. A fixed bin width of 25 HU is used for grayscale discretization, and a symmetric grayscale co-occurrence matrix is ​​enabled. This yields the normalized voxel matrix V and the binary mask M registered with it. The above values ​​constitute a set of reproducible specific parameter combinations. In other imaging modalities or training data, the isotropic voxel spacing, grayscale upper and lower limits, and discretization parameters can be determined according to the voxel spacing distribution and grayscale distribution of the training samples and remain unchanged during the application stage.

[0030] Regions of interest (ROI) masks can be obtained through manual delineation, semi-automatic region growing, or trained segmentation networks. As a quality control method, two observers can each delineate portions of the samples, and samples with boundary differences exceeding 5 mm and appearing consecutively in three or more layers are submitted to a third observer for review. During the application phase, only the generated mask files are read; operations performed on the human body are not included. Example 3: Multi-path radiomics feature extraction

[0031] like Figure 3 As shown, the normalized voxel matrix V and the mask M are input into the feature extraction module 40. The original image branch 41 directly calculates first-order statistics, three-dimensional shape, and texture features on the voxel set defined by M; the spatial filtering branch 42 uses a three-dimensional Gaussian Laplacian filter with a scale parameter σ = 3 mm to obtain the filtered volume data, and then calculates first-order statistics and texture features; the multi-scale transformation branch 43 uses three-dimensional wavelet decomposition to obtain LLL, LLH, LHL, LHH, HLL, HLH, HHL, and HHH subbands, and extracts first-order statistics and texture features from each subband.

[0032] In a reproducible implementation, using PyRadiomics 3.0.1-compatible feature definitions, a total of 957 candidate features were generated, including 120 original image features, 93 Gaussian Laplacian filter features, and 744 wavelet transform features. Texture features may include gray-level co-occurrence matrix, gray-level run-length matrix, gray-level size region matrix, adjacent gray-level tonal difference matrix, and gray-level correlation matrix features. The system records the name, filter branch, parameters, mask identifier, and calculation version for each feature for verification.

[0033] For other implementations, multiple Gaussian scale parameters can be set, or spatial transformations such as local binary mode, Gabor filtering, and discrete cosine transform can be used, as long as the output is still organized according to a fixed feature dictionary and consistent transformation parameters are used in the training and application phases. Example 4: Cascaded Filtering and Model Training

[0034] like Figure 4 As shown, model training is performed in a training environment isolated from the application phase. The example dataset includes 164 de-labeled samples, with 100 first-label samples and 64 second-label samples. A fixed random seed of 42 is used for stratified random partitioning, resulting in a training set of 114 samples (70 first-label samples) and a validation set of 50 samples (30 first-label samples). The labels are generated from historical data and are used solely for offline supervised learning.

[0035] First, the stability screening unit 51 calculates the intra-group correlation coefficients of candidate features from repeated segmentation samples or multi-observer segmentation samples, and retains features based on a stability threshold determined within the training samples. In this specific implementation, intra-group correlation coefficients are calculated from 30 repeated segmentation samples, and features with intra-group correlation coefficients greater than 0.75 are retained. Second, only the mean μj and standard deviation σj of each feature are fitted using the training set, and a Z-score transformation is performed according to zj = (xj - μj) / σj, while removing features with near-zero variance. The trained μj and σj are fixed in the model configuration file, and the validation set and application data directly call this configuration without refitting.

[0036] Then, the sparsity screening unit 52 determines the regularization parameters of the regression model with regularization terms through K-fold cross-validation within the training samples and retains features with non-zero coefficients; the correlation screening unit 53 traverses candidate features from high to low importance and removes redundant features based on the correlation threshold determined within the training samples. In this embodiment, the minimum absolute shrinkage and selection operator of ten-fold cross-validation is used for sparsity screening, and the regularization parameters are selected according to a standard error rule. The 15 candidate features with higher importance among the non-zero coefficients are retained first; for new features, if the absolute value of the Spearman correlation coefficient between it and any of the already retained features is not less than 0.75, it is marked as redundant, otherwise it is retained. Finally, the total number of variables in the fusion model is controlled according to the constraint that the ratio of the number of positive events to the number of predictive variables is not less than 10; when at least one structured attribute variable is included, the number of target image features is limited to no more than six.

[0037] The target feature subspace may include at least two target image features determined from the candidate feature set by the cascaded filtering results of training samples. The six target image features obtained in this embodiment are, in order: f1, log-sigma-3-0-mm-3D_glszm_SmallAreaEmphasis; f2, log-sigma-3-0-mm-3D_glszm_GrayLevelNonUniformity; f3, original_shape_Maximum2DDiameterRow; f4, log-sigma-3-0-mm-3D_firstorder_InterquartileRange; f5, diagnostics_Image-original_Mean; f6, wavelet-LLL_glcm_InverseVariance. Here, LLL represents the wavelet subband obtained by applying low-pass filtering in all three directions of the three-dimensional wavelet decomposition. The median intra-group correlation coefficient of the six features is 0.857, ranging from 0.844 to 0.883.

[0038] For an implementation containing m target image features, the radiomics score can be calculated as R = β0 + β1 × f1 + ... + βm × fm, where m is a positive integer not less than 2, and β0 to βm are fitted and fixed only by training samples. In this embodiment, m = 6, and the score is calculated as follows: Rad-score = 0.6140 + 0.0718 × f1 - 0.0475 × f2 + 0.1986 × f3 + 0.1639 × f4 + 0.1440 × f5 + 0.1995 × f6. The feature values ​​in the formula are all values ​​after Z-score transformation using the training set parameters. If the feature implementation library names f5 as original_firstorder_Mean, its mapping relationship with the training field should be saved in the model configuration file and used after passing the numerical consistency test. The specific coefficients are used to provide a reproducible narrow-layer implementation; different weight parameters can be obtained for other training data.

[0039] The radiomics score or target image feature vector is concatenated with age coding and a pre-determined image infiltration level coding. Structured attributes are read only as existing data fields; this invention does not generate the image infiltration level. As a simplified linear branch, age, image infiltration level, and Rad-score can be used to train a logistic regression model; the example regression coefficients are -0.04, 1.29, and 7.45, respectively, and the intercept is obtained through maximum likelihood estimation of the training set. These regression coefficients are also specific model parameters determined by the corresponding training samples.

[0040] As a relatively complete fusion model, six target image features, age coding, and image infiltration level coding can be used as inputs to train logistic regression, support vector machine, random forest, extreme gradient boosting, and lightweight gradient boosting models, respectively. Each model employs class-balanced weights and performs 10-fold cross-validation on the training set. During validation, the normalization parameters, feature set, or model parameters are not adjusted. The model file must at least save the feature order, missing value rules, normalization parameters, model weights, training timestamps, and version checksums. Example 5: Dual-path classification and interpretable data packets

[0041] like Figure 5 As shown, the same fused feature vector is input into the linear classification branch 61 and the nonlinear ensemble classification branch 62, respectively. In this embodiment, logistic regression is used as the linear classification branch 61, and random forest is used as the nonlinear ensemble classification branch 62. Each branch outputs a classification reference value between 0 and 1. One branch can be selected as the primary numerical output based on pre-defined validation rules, while the value of the other branch is retained for difference indication; alternatively, in other implementations supported by the specification, validated weights can be used for weighting.

[0042] The contribution calculation unit 71 uses the training set as background samples to calculate the Shapley Additive Explanations (SHAP) values ​​for the two branches. For a sample x, the model output is represented as g(x) = φ0 + Σφj, where φ0 is the expected background output and φj is the contribution value of the j-th feature component. The global contribution is represented by the mean of |φj| over the sample set, and the local contribution vector is composed of the φj values ​​of the current sample.

[0043] The interaction strength calculation unit 72 calculates the bivariate partial dependency function for the nonlinear ensemble classification branch 62. For features j and k, the difference between the centered bivariate partial dependency and the two univariate partial dependencies can be taken as the non-additive interaction term, and the Friedman H statistic is obtained by taking the square root of the ratio of the sum of squares of this interaction term to the sum of squares of the bivariate partial dependencies. The computation grid resolution can be set to 50, and the trapezoidal rule can be used for numerical integration. The closer the H statistic is to 0, the weaker the non-additive interaction in the model.

[0044] The data packet generation unit 73 writes the main classification reference values, the difference between the reference values ​​of the two branches, the global contribution ranking, the local contribution vector of the current sample, the Friedman H statistic matrix, preprocessing parameters, model version, and timestamp into a structured data packet. For easy machine reading, the data packet can be in JSON, XML, DICOM structured report, or database record format; for easy manual review, contribution bar charts, decision path diagrams, and interaction matrix diagrams can be generated simultaneously.

[0045] In one optional implementation, the model difference index D is defined as |pL - pN| + λ × [1 - cos(ΦL, ΦN)] / 2, where pL and pN are the classification reference values ​​for the linear and nonlinear branches, respectively, ΦL and ΦN are the normalized local contribution vectors, and λ is the non-negative weight determined on the validation set. When D exceeds a preset threshold, the output module 80 displays a technical prompt "The model difference is large, and it is recommended to review the data quality," without outputting disease diagnosis or treatment suggestions. Example 6: Electronic Device and Program Implementation

[0046] like Figure 6 As shown, the electronic device 90 includes a processor 91, a memory 92, a medical image data interface 93, a display interface 94, a network interface 95, and a bus 96. The processor 91 may be a central processing unit, a graphics processor, a neural network processor, or a combination thereof; the memory 92 may include volatile memory and non-volatile memory; the medical image data interface 93 may be connected to an image archiving and communication system.

[0047] The memory 92 stores the preprocessing configuration, feature dictionary, normalization parameters, classification model, and interpreter. After acquiring data from the medical image data interface 93, the processor 91 sequentially performs image standardization, feature extraction, target feature mapping, fusion classification, and interpretation calculation. For multiple filtering branches, the processor 91 can establish a parallel task queue; for wavelet volume data, the texture matrix can be calculated in sub-band blocks to reduce peak storage usage.

[0048] The following exception control measures can be implemented during program execution: When the spatial matrices of the image or mask are inconsistent, feature extraction is stopped and the error is recorded; when target features are missing, they are filled in using the median or other preset values ​​fixed during training, and the filling identifier is recorded in the data packet; when the order of input fields or model validation values ​​are inconsistent, classification mapping is not performed. These controls help reduce untraceable outputs caused by inconsistent inputs.

[0049] The media hosting the above programs can be solid-state drives, hard disks, read-only memory, flash memory, optical discs, or other non-transitory computer-readable storage media. The programs can also be deployed on access-controlled local servers, dedicated workstations, or compliant cloud computing environments. Example 7: Verification Example

[0050] To verify the above processing flow, the area under the receiver operating characteristic (ROC) curve, calibration curve, Brier score, precision-recall curve, and mean precision were calculated on both the training and validation sets. Confidence intervals were estimated using 1000 bootstrapping resampling operations. The DeLong test can be used for model comparison, but this statistical analysis is only used to verify the algorithm output and does not constitute a step in the method of this invention to diagnose an individual.

[0051] The area under the curve (AUC) for the model using only structured attributes on the training and validation sets is approximately 0.753 and 0.794, respectively; for the radiomics model using six target image features, it is approximately 0.919 and 0.813, respectively; for the fused logistic regression model, it is approximately 0.907 and 0.868, respectively; for the fused random forest model, it is approximately 0.839 and 0.857, respectively; and for the simplified model, it is approximately 0.846. The difference in AUC between the fused logistic regression model and the structured attribute model is approximately 0.086, corresponding to a p-value of approximately 0.048 in this validation example.

[0052] At a preset threshold, the sensitivity of the fusion model on the validation set is approximately 0.833, and the specificity ranges from approximately 0.750 to 0.800 depending on the threshold setting. These values ​​are intended to illustrate that the feature processing and model fusion can produce discriminative computer output, and do not necessarily represent the same results obtained on other datasets, devices, or populations.

[0053] SHAP analysis showed that `original_shape_Maximum2DDiameterRow` and `wavelet-LLL_glcm_InverseVariance` had high mean absolute contributions in both classification branches. Exemplary interactions detected in the random forest branches included: H = 0.136 between `original_shape_Maximum2DDiameterRow` and `log-sigma-3-0-mm-3D_glszm_SmallAreaEmphasis`; H = 0.109 between `original_shape_Maximum2DDiameterRow` and image infiltration level; H = 0.101 between `original_shape_Maximum2DDiameterRow` and age; and H = 0.100 between age and image infiltration level. This demonstrates that the interaction matrix can provide model behavior information that is not available when only a single classification value is output. Example 8: Replaceable Implementation

[0054] This invention is not limited to gastric CT images. While maintaining the overall concept of pre-acquired 3D digital images, registration masks, multi-path feature extraction, cascaded filtering, fusion classification, and interpretation output, the input can be replaced with CT, MRI, positron emission tomography (PET) images of other body parts, or combinations thereof. For different imaging modalities, the grayscale limiting range, resampling interval, and discretization parameters can be determined by the training set and fixed along with the model configuration.

[0055] Linear classification branches can be implemented using linear discriminant models, generalized linear models, or linear support vector machines; nonlinear classification branches can be implemented using random forests, gradient boosting trees, kernel support vector machines, or neural networks. Contribution algorithms can also employ other algorithms that satisfy local additivity or counterfactual interpretation requirements, while interaction algorithms can use permutation interactions, SHAP interaction values, or variance decomposition-based indices.

[0056] It should be understood that the single numerical values, thresholds, number of features, and model coefficients described in this specification are used to illustrate a reproducible specific implementation. Under the technical concept of maintaining parameters determined during the training phase and fixed with each model version, consistent invocation during the application phase, and multi-path extraction, cascaded filtering, fusion classification, and interpretation output, those skilled in the art can make adaptive adjustments to module merging, execution order, parameter selection, and model implementation. The specific scope of protection is determined by the claims.

Claims

1. A medical image processing method integrating radiomics and machine learning, characterized in that, The method, executed by a computer device on pre-acquired and stored medical image data, includes: Acquire the three-dimensional medical image data to be processed and the region of interest mask spatially registered with the three-dimensional medical image data; Voxel resampling, grayscale limiting, and grayscale discretization are performed on the three-dimensional medical image data to obtain a standardized voxel matrix; Using the region defined by the region of interest mask as the computation domain, a parallel feature extraction path is constructed, including the original image branch, the spatial filtering branch, and the multi-scale transformation branch. Statistical features, shape features, and texture features are extracted from the parallel feature extraction path to form an image omics feature vector. The radiomics feature vector is normalized based on normalization parameters determined only by training samples, and the normalized radiomics feature vector is mapped to the target feature subspace in sequence according to stability constraints, sparsity constraints, and correlation constraints. The target image features in the target feature subspace are encoded and concatenated with the pre-acquired structured attribute features to obtain a fused feature vector. The fused feature vector is then input into a trained classification model to obtain a classification reference value that characterizes the degree of matching between the three-dimensional medical image data and the predetermined image category. Based on the output of the classification model, the contribution of each feature component and the interaction strength of at least one pair of feature components are calculated, and an interpretable data packet including the classification reference value, contribution ranking and interaction data is generated. The interpretable data packet is then output to a storage device or a display device.

2. The medical image processing method according to claim 1, characterized in that, The three-dimensional medical imaging data includes computed tomography volume data, magnetic resonance imaging volume data, positron emission tomography volume data, or a spatial registration combination of the aforementioned volume data. The voxel resampling includes resampling the three-dimensional medical image data to a preset isotropic voxel spacing; the grayscale limiting includes limiting the voxel value between a preset lower grayscale limit and a preset upper grayscale limit; and the grayscale discretization adopts a fixed bin width method or a fixed number of bins method. The isotropic voxel spacing, preset lower gray level, preset upper gray level, and gray level discretization parameters are predetermined based on the imaging modality, voxel spacing distribution, and gray level distribution of the training samples, and are stored in association with the model version so that they can be called consistently during the application phase.

3. The medical image processing method according to claim 1, characterized in that, The original image branch is used to extract first-order statistical features and shape features from the unfiltered normalized voxel matrix; The spatial filtering branch is used to perform three-dimensional Gaussian Laplacian filtering on the normalized voxel matrix and extract first-order statistical features and texture features from the filtering results. The multi-scale transformation branch is used to perform three-dimensional wavelet decomposition on the normalized voxel matrix and extract texture features from multiple wavelet subbands; The texture features include at least one of the following: gray-level co-occurrence matrix features, gray-level run-length matrix features, gray-level size region matrix features, adjacent gray-level tone difference matrix features, and gray-level correlation matrix features.

4. The medical image processing method according to claim 1, characterized in that, The normalization parameters are obtained solely from fitting the training samples; The stability constraint includes calculating the intra-group correlation coefficient based on repeatedly segmented samples or multi-observer segmented samples, and retaining the feature that the intra-group correlation coefficient reaches a preset stability threshold. The sparsity constraint includes removing near-zero variance features and using a regression model with a regularization term to filter non-zero coefficient features, wherein the regularization parameter of the regression model with a regularization term is determined by performing cross-validation within the training samples. The correlation constraint includes traversing the features retained by the sparsity constraint in order of feature importance, and removing features whose absolute value of the correlation coefficient with the retained features reaches a preset correlation threshold. The preset stability threshold, regularization parameter, and preset correlation threshold are all determined within the training samples, and the dimension of the target feature subspace is limited according to the ratio of the number of positive events to the number of candidate predictive variables in the training samples.

5. The medical image processing method according to claim 4, characterized in that, The target feature subspace includes m target image features determined from the candidate feature set by the stability constraint, sparsity constraint and correlation constraint, where m is a positive integer not less than 2; The candidate feature set includes: small region emphasis features of gray-scale region matrix after three-dimensional Laplace Gaussian filtering, gray-scale non-uniformity features of gray-scale region matrix after three-dimensional Laplace Gaussian filtering, maximum two-dimensional path features of the original image, first-order interquartile range features after three-dimensional Laplace Gaussian filtering, mean features of the original image, and inverse variance features of LLL subband gray-scale co-occurrence matrix obtained by three-dimensional wavelet transform with low-pass filtering in all three spatial directions; Based on the normalized m target image features, a radiomics score R is constructed, R = β0 + β1 × f1 + ... + βm × fm, where f1 to fm are the m target image features, and β0 to βm are weight parameters that are fitted only based on training samples and are fixed with the model version.

6. The medical image processing method according to claim 5, characterized in that, m=6, the six target image features are denoted as f1 to f6 in the order listed in claim 5, and the image omics score is: Rad-score = 0.6140 + 0.0718×f1 - 0.0475×f2 + 0.1986×f3 + 0.1639×f4 + 0.1440×f5 + 0.1995×f6.

7. The medical image processing method according to claim 1, characterized in that, The structured attribute features include age coding and pre-determined image infiltration level coding; The classification model includes a linear classification branch and a nonlinear ensemble classification branch. The linear classification branch includes a logistic regression model or a support vector machine model. The nonlinear ensemble classification branch includes a random forest model, an extreme gradient boosting model, or a lightweight gradient boosting model. Each classification branch is trained using class-balanced weights, and the main output branch is determined based on at least two of the following metrics: discrimination metric, calibration metric, and precision-recall metric obtained from cross-validation. The contribution is obtained by calculating the Shapley additive explanatory value of the two classification branches based on the same background sample set, and the contribution ranking is determined according to the mean of the absolute values ​​of the Shapley additive explanatory values ​​of each feature component. The interaction strength is obtained by calculating the Friedman H statistic for the bivariate partial dependency function of the nonlinear ensemble classification branch; The interpretable data package also includes the difference in classification reference values ​​between the two classification branches, the local contribution vector, and the Friedman H statistic matrix.

8. A medical image processing system integrating radiomics and machine learning, characterized in that, include: The data receiving module is used to read pre-acquired and stored 3D medical image data and region of interest masks; The image preprocessing module is used to perform voxel resampling, grayscale limiting, and grayscale discretization on the three-dimensional medical image data based on the preprocessing parameters stored in association with the model version. The feature extraction module is used to extract radiomics feature vectors through the original image branch, spatial filtering branch, and multi-scale transformation branch; The feature selection module is used to generate target image features based on normalization parameters, stability constraints, sparsity constraints, and correlation constraints determined by training samples. The fusion classification module is used to encode and concatenate the target image features and structured attribute features, and then input them into the trained classification model to generate classification reference values; The interpretation output module is used to calculate the contribution and interaction strength of the feature components, generate and output an interpretable data package including the classification reference value, contribution ranking, interaction data, model version and preprocessing parameters.

9. An electronic device, comprising a processor, a memory, a medical imaging data interface, and a display interface, wherein the processor is communicatively connected to the memory, the medical imaging data interface, and the display interface, characterized in that, The memory stores a computer program, which, when executed by the processor, causes the electronic device to perform the medical image processing method according to any one of claims 1 to 7.

10. A non-transitory computer-readable storage medium having a computer program stored thereon, characterized in that, When the computer program is executed by a processor, it implements the medical image processing method according to any one of claims 1 to 7.

Images