Tumor area image texture analysis method and system

Abstract

The invention discloses a tumor area image texture analysis method and system, and relates to the technical field of image processing. Comprising the following steps: S1, image preprocessing: preprocessing obtained multi-modal data to obtain preprocessed multi-modal data and a feature image; s2, biological association: calibrating the feature image, a multi-center equipment parameter library and reference standard data, and obtaining the association between a biological marker and a texture parameter through a machine learning model; s3, prediction modeling: obtaining a risk score according to the unified image domain and the real-time monitoring data; and S4, generating a treatment scheme: according to the risk score and the real-time surgical navigation data, adjusting a treatment period through reinforcement learning, and obtaining a personalized treatment scheme. According to the method, the problems of micro-feature loss, data deviation and poor clinical interpretability in traditional texture analysis are solved.

Description

Technical Field

[0001] The present invention relates to the field of image processing technology, and in particular to a method and system for analyzing texture of tumor region images. Background Art

[0002] With the continuous advancement of medical imaging technology, medical images have become an important basis for the clinical diagnosis of tumor lesions. Modalities such as computed tomography (CT), magnetic resonance imaging (MRI), and positron emission tomography (PET) are widely used for tumor detection, staging, and therapeutic efficacy assessment. However, traditional imaging diagnosis relies primarily on the subjective judgment of radiologists, which has certain limitations, such as low reading efficiency and large individual variability.

[0003] In recent years, the integration of medical image processing and artificial intelligence technologies has driven the development of quantitative image analysis. Image texture analysis, as a key method for extracting image features, has shown great potential in tumor region identification and classification. Texture information can reflect the spatial distribution of internal tissue structures, helping to reveal subtle changes in lesions that are difficult to detect with the naked eye. This provides an objective basis for early tumor detection, determination of benign or malignant status, and prognostic prediction.

[0004] Currently, a variety of texture analysis methods have been applied to tumor imaging research, including gray-level co-occurrence matrix (GLCM), gray-level run-length matrix (GLRLM), gray-level gradient co-occurrence matrix (GLGSM), wavelet transform, local binary pattern (LBP), and deep learning-based texture feature extraction methods. These methods quantify the spatial relationship and intensity changes between pixels or voxels in an image to construct discriminative feature vectors, which are then used for training and classification of machine learning or deep learning models.

[0005] The existing Chinese invention patent with publication number CN112768072A discloses a cancer clinical indicator evaluation system based on a qualitative imaging genomics algorithm, which includes: a data input module for inputting multi-dimensional patient clinical data, wherein the multi-dimensional patient clinical data includes: the patient's preoperative multimodal medical imaging atlas and clinical data; a feature extraction and optimization module, connected to the data input module, for extracting imaging genomics quantitative feature values ​​and clinical data features, and performing dimensionality reduction and optimization on the quantitative feature values; a qualitative model construction module, connected to the feature extraction and optimization module, for constructing an evaluation model for predicting various clinical indicators of patient cancer through a qualitative imaging genomics algorithm; an evaluation module, connected to the qualitative model construction module, for providing non-invasive, accurate and non-subjective evaluation of early diagnosis, prognosis and drug efficacy subtypes of cancer.

[0006] However, existing dynamic texture analysis methods employ smoothing constraints to maintain image registration stability during submillimeter precision recognition and clinical interpretability, which inevitably filters out microtexture features with early diagnostic value (such as the directional alignment pattern of collagen fibers <0.5 mm after immunotherapy). Furthermore, due to the lack of a standardized multicenter biomarker validation framework, the correlation between dynamic texture parameters (such as the spatiotemporal heterogeneity index) and key biological processes in the tumor microenvironment (such as metabolic remodeling in hypoxic regions) is difficult to reliably verify. This can lead to systematic deviations in feature parameters acquired using different imaging devices and scanning protocols, which in turn increases the rate of misjudgment of treatment responses in clinical practice and hinders the large-scale clinical application of personalized cancer diagnosis and treatment plans. Summary of the Invention

[0007] The purpose of the present invention is to provide a method and system for analyzing texture of tumor region images to solve the problems raised in the above background technology.

[0008] To achieve the above objectives, the present invention provides the following technical solutions: 1. A method for analyzing texture of a tumor region image, comprising:

[0009] S1: Image preprocessing: preprocessing the acquired multimodal data to obtain preprocessed multimodal data, and performing feature extraction on the preprocessed multimodal data to obtain a feature image;

[0010] S2: Biological Correlation: Calibrate the feature image with a multi-center device parameter library and reference standard data, and obtain the correlation between biological markers and texture parameters through a machine learning model, including:

[0011] S2.1: Data preprocessing: Calibrate the constructed multi-center equipment parameter library and reference standard data to obtain a unified image domain;

[0012] S2.2: Biological Validation: Combining the unified image domain with the biomarker monitoring data to obtain a fused data matrix, and using the fused data matrix as input to a random forest model to output correlations between texture parameters and biomarkers;

[0013] S3: Predictive modeling: Based on the unified image domain and real-time monitoring data, a spatiotemporal heterogeneity-immune clonal diversity joint model is constructed to obtain a risk score;

[0014] S4: Generate treatment plan: Based on the risk score and real-time surgical navigation data, adjust the treatment cycle through reinforcement learning to obtain a personalized treatment plan.

[0015] Furthermore, the feature image is obtained, including:

[0016] S1.1: Image processing: Align the T1-weighted image, elastic modulus map, and CD8+ immunohistochemical slice image to the same coordinate system using predefined anatomical landmarks. Resample the voxel sizes of the elastic modulus map and CD8+ immunohistochemical slice image to match those of the T1-weighted image. Obtain a binary mask matrix for the elastic modulus map and an immune heat map for the CD8+ immunohistochemical slice image.

[0017] S1.2: Microtexture extraction: Based on the T1-weighted image, binary mask matrix and immune heat map, a three-dimensional multimodal data cube is constructed, and the three-dimensional multimodal data cube is processed using a gray-level co-occurrence matrix, local binary pattern and Gabor filter to obtain a multimodal feature matrix and a three-dimensional fusion image.

[0018] Furthermore, the three-dimensional fusion image includes a background color and an overlay color. The background color is set to the grayscale of the T1 weighted image, and the overlay color is set to different colors according to the LBP peak area, Gabor high response area and CD8+ high density area.

[0019] Furthermore, obtaining the binary mask matrix of the elastic modulus map and the immune heat map of the CD8+ immunohistochemistry slice image includes:

[0020] S1.1.1: Data normalization: Based on the T1-weighted images, register the elastic modulus map and CD8+ immunohistochemistry slice images to the same coordinate system using affine transformation, and downsample the CD8+ immunohistochemistry slice images;

[0021] S1.1.2: Generate elastic modulus mask: Compare the magnitude of each modulus in the standardized elastic modulus map with a preset modulus threshold, and mark each spatial position in the standardized elastic modulus map based on the comparison result, specifically:

[0022] When the modulus at the spatial position is smaller than the preset modulus threshold, the spatial position is marked as 1; otherwise, the spatial position is marked as 0;

[0023] S1.1.3: CD8+ immune markers: Based on the CD8+ cell localization coordinates and the global coordinate system of the T1-weighted image, obtain the nuclear density size, compare the nuclear density size with the preset cell density threshold, and delete the center coordinates of the immune cells based on the comparison results, specifically:

[0024] When the nuclear density is smaller than the preset cell density threshold, the central coordinates of the immune cells corresponding to the nuclear density are deleted; otherwise, the central coordinates of the immune cells corresponding to the nuclear density are retained.

[0025] Furthermore, a unified image domain is obtained, including:

[0026] S2.1.1: Parameter library standardization: Obtain standardized characteristic values ​​of the multi-center device parameter library based on the data values ​​in the multi-center device parameter library and the characteristic means and standard deviations of the current device / reference device in healthy tissue;

[0027] S2.1.2: Protocol Conversion: Perform structural modeling based on the liver module and tumor classification ratio to obtain a customized digital phantom. Simultaneously, obtain the calibration slope and intercept using the sequence parameters of Center A and Center B. Based on the calibration slope and intercept, determine the standardized signal value, specifically:

[0028]

[0029] in: is the normalized signal value, is the original device signal value, is the calibration slope, is the calibration intercept;

[0030] S2.1.3: Image domain unification: The standardized signal values ​​and the dual-center paired DICOM files are used as inputs to the CycleGAN model, and the synthesized feature images are obtained as output.

[0031] Furthermore, the outputs obtained the correlation between texture parameters and biological markers, including:

[0032] S2.2.1: Spatial transcriptome matching: Based on the DAPI-stained image, align the transcriptome array with the pathology slide, and obtain the mean and standard deviation of the image area covered by the transcriptome spatial unit to determine the spatial variation characteristics of the tumor microenvironment. Specifically:

[0033]

[0034] in: is the heterogeneity index, is the characteristic measurement value of the m-th spatial unit, is the characteristic average of all spatial units, is the global standard deviation, is the total number of spatial units analyzed, is the spatial unit index;

[0035] S2.2.2: CTC proteomic association analysis: Spearman correlation analysis was used to temporally correlate CTC characteristics with imaging markers to obtain the Spearman rank correlation coefficient.

[0036] Furthermore, a risk score is obtained, including:

[0037] S3.1: Data Fusion: Each patient is set as a central node, and the connection weights between adjacent nodes are set based on the apparent diffusion coefficient heterogeneity index of the unified image domain. At the same time, the set central nodes are connected through the GNN model, and the feature vectors of each node at each layer are determined based on the connection weights;

[0038] S3.2: Constructing a hierarchical model: Based on the treatment response probability of the GNN model and the actual treatment response, the target loss function is obtained, specifically:

[0039]

[0040] in: is the target loss function, is the true label of the rth sample, is the model-predicted probability of treatment response for the rth sample, is the total number of samples, is the sample index, is the tree structure complexity penalty coefficient, is the total number of leaf nodes, is the L2 regularization coefficient, is the L2 norm squared of the weight.

[0041] Go further and get a personalized treatment plan, including:

[0042] S4.1: Toxicity prediction: Based on the elastic modulus graph and ADC graph, obtain the elastic modulus standard deviation, and determine the probability of radiation pneumonitis based on the elastic modulus standard deviation. At the same time, classify the toxicity level based on the probability of radiation pneumonitis;

[0043] S4.2: Solution optimization: Based on the risk score and toxicity level, a reward function is obtained through a reinforcement learning strategy, and the current solution is optimized based on the reward function and the Q-learning algorithm.

[0044] Furthermore, the probability of occurrence of radiation pneumonitis is compared with a preset probability threshold range, and a corresponding clinical response is determined based on the comparison result, specifically:

[0045] When the probability of radiation pneumonitis is less than the lower limit of the preset probability threshold, the current radiotherapy plan is continued; when the probability of radiation pneumonitis is within the preset probability threshold, the frequency of high-resolution CT monitoring is increased; when the probability of radiation pneumonitis is greater than the upper limit of the preset probability threshold, the dose is reduced and serum testing is performed.

[0046] A tumor region image texture analysis system uses any one of the above-mentioned tumor region image texture analysis methods.

[0047] Compared with the prior art, the present invention has the following beneficial effects:

[0048] First, the present invention fuses T1-weighted images, elastic modulus maps, and CD8+ immune heatmaps, and extracts microtexture features using gray-level co-occurrence matrices, local binary patterns, and Gabor filters, enabling submillimeter tumor feature recognition. Furthermore, the Spearman rank correlation coefficient is used to analyze the temporal association between circulating tumor cell counts and imaging markers, ensuring a strong correlation between texture parameters and clinical markers.

[0049] Second, the present invention standardizes the device parameter library by using the mean and standard deviation of healthy tissue features, thereby eliminating the systematic deviation between different devices. Furthermore, the CycleGAN model is used to synthesize the dual-center DICOM files into standardized feature images, thereby resolving the problem of protocol differences between devices.

[0050] Third: The present invention uses patients as nodes to construct a GNN hierarchical model, quantifies spatiotemporal heterogeneity through connection weights, and predicts the probability of treatment response in combination with the target loss function. At the same time, the corresponding probability of radiation pneumonia is obtained through the standard deviation of the elastic modulus and the ADC coefficient of variation, and the toxicity level is divided to optimize the treatment plan. This can not only intuitively display the tumor boundary characteristics and assist surgical navigation, but also dynamically adjust the treatment cycle to balance efficacy and safety. BRIEF DESCRIPTION OF THE DRAWINGS

[0051] Figure 1 Schematic diagram of the process of the tumor region image texture analysis method of the present invention;

[0052] Figure 2 This is a diagram of multimodal data processing in the present invention;

[0053] Figure 3 It is the LBP mode distribution diagram of the present invention;

[0054] Figure 4 It is the T1-weighted image in the present invention;

[0055] Figure 5 is the elastic modulus diagram of the present invention;

[0056] Figure 6 This is the three-dimensional fusion image in the present invention. DETAILED DESCRIPTION

[0057] The following will clearly and completely describe the technical solutions in the embodiments of the present invention in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by ordinary technicians in this field without making creative efforts are within the scope of protection of the present invention.

[0058] Existing dynamic texture analysis methods use smoothing constraints to maintain image registration stability during submillimeter precision recognition and clinical interpretability, which inevitably filters out microtexture features with early diagnostic value (such as collagen fiber orientation patterns <0.5mm after immunotherapy). Furthermore, due to the lack of a standardized multicenter biomarker validation framework, the correlation between dynamic texture parameters (such as spatiotemporal heterogeneity index) and key biological processes in the tumor microenvironment (such as metabolic remodeling in hypoxic areas) is difficult to reliably verify. This can lead to systematic deviations in feature parameters obtained under different imaging devices and scanning protocols, which not only increases the misjudgment rate of treatment response in clinical practice but also restricts the large-scale clinical application of personalized cancer diagnosis and treatment plans. The technical solution of this application constructs a three-dimensional multimodal data cube using T1-weighted images, elastic modulus maps, and CD8+ immunohistochemistry slice images, generating a multimodal feature matrix and a three-dimensional fusion map. It also utilizes a multicenter device parameter library and reference standard data to unify the image domain, eliminate device differences, and establish the correlation between texture parameters and biomarkers using spatial transcriptome and CTC proteomics data. At the same time, based on the unified image domain and real-time monitoring data, the risk of tumor microenvironment heterogeneity is quantified, and the treatment cycle and dose are dynamically adjusted to generate personalized treatment plans, thus solving the problems of micro-feature loss, data bias and poor clinical interpretability in traditional texture analysis.

[0059] Example 1

[0060] refer to Figures 1-6 This embodiment provides a method for analyzing the texture of a tumor region image. The method specifically includes the following steps:

[0061] Step S1: Image preprocessing. The obtained multimodal data is preprocessed to obtain the preprocessed multimodal data, and the preprocessed multimodal data is subjected to feature extraction to obtain the corresponding feature image. The details are as follows:

[0062] Step S1.1: Image processing. This involves acquiring T1-weighted images, elastic modulus maps, and CD8+ immunohistochemical slice images using the MRI equipment, magnetic resonance elastography system, and pathology scanning system. These images are then aligned to the same coordinate system using pre-defined anatomical landmarks. Furthermore, the voxel sizes of the elastic modulus maps and CD8+ immunohistochemical slice images are resampled to match those of the T1-weighted images.

[0063] Furthermore, in the elastic modulus map, each spatial location in the elastic modulus map is marked according to a preset modulus threshold, resulting in a binary mask matrix of the same size as the original image. Simultaneously, the CD8+ cell count in the CD8+ immunohistochemistry slice image is converted to standard units, and a corresponding set of coordinate points is obtained according to a preset cell density threshold. Based on this set of coordinate points, a corresponding immune heat map is generated.

[0064] Step S1.2: Microtexture Extraction. Based on the T1-weighted images, elastic modulus maps, and immune thermograms obtained in step S1.1, a three-dimensional multimodal data cube is constructed. Each voxel in the three-dimensional multimodal data cube includes grayscale value, elastic modulus, and CD8+ cell density. Furthermore, the constructed three-dimensional multimodal data cube undergoes non-local mean denoising using a preset search window, similarity block, and filter strength. Simultaneously, balanced denoising and edge preservation are employed, with an iterative process of 0.05 time steps. Specifically, smoothing is enhanced in uniform regions with elastic modulus gradients less than 3 kPa / mm, while diffusion is suppressed in high-gradient regions at the tumor boundary. Furthermore, the CT data are grayscale normalized with a window width of 400 HU and a window level of 40 HU. The MRI data are corrected using an N4 bias field, histogram matched, and normalized to the standard brain atlas grayscale range, thereby obtaining the preprocessed three-dimensional multimodal data cube.

[0065] In this embodiment, the pre-processed three-dimensional multimodal data cube is processed by gray-level co-occurrence matrix, local binary pattern and Gabor filter, and a multimodal feature matrix is ​​obtained, such as Figure 2 Specifically, the modulus values ​​in the elastic modulus diagram are discretized at intervals of 5 kPa, and the corresponding energy, contrast, correlation, variance, inverse gap, and entropy data are obtained based on the discretized modulus values.

[0066] Furthermore, the original T1 weighted image is processed by LBP and Gabor filtering. Specifically, the radius is set to 2 pixels, the 8-neighborhood uniform pattern is used, and the processing is performed in a rotation-invariant manner, and the corresponding LBP pattern distribution map is output, as shown in FIG. Figure 3As shown. Furthermore, the constructed Gabor filter is used to extract the modulus mean, standard deviation, and energy of each scale / direction combination in the original T1-weighted image. In this embodiment, the constructed Gabor filter includes three scales, and the number of directions of each scale is set to 8. The bandwidth of adjacent scales increases by two times, that is, the bandwidth of the first scale is set to 0.8, the bandwidth of the second scale is set to 1.6, and the bandwidth of the third scale is set to 3.2. At the same time, the wavelengths of adjacent scales increase in sequence, that is, the wavelength of the first scale is set to 0.5 mm, the wavelength of the second scale is set to 1 mm, and the wavelength of the third scale is set to 2 mm.

[0067] Furthermore, a deep learning model (such as Mask R-CNN) is used to label CD8+ cell nuclei in CD8+ immunohistochemistry slice images and obtain the number of cells within each grid. Specifically, the tumor ROI is divided into multiple grids, and the corresponding density ratio of CD8+ cells in each grid is determined. Based on the density ratio, the energy characteristics of the corresponding CD8+ cell distribution are obtained, specifically:

[0068]

[0069] in: is the energy characteristic of CD8+ cell distribution, is the density ratio of the jth partition, is the CD8+ cell density of the jth partition, is the total number of spatial partitions, is the current partition index, is the total density of CD8+ cells in all partitions.

[0070] In this embodiment, after processing the three-dimensional multimodal data cube using a gray-level co-occurrence matrix, local binary pattern, Gabor filter, and energy feature extraction, each feature data obtained is standardized, and based on the standardized feature data, a corresponding multimodal feature matrix and a corresponding three-dimensional fusion image are constructed. Specifically, in the obtained three-dimensional fusion image, the T1 image grayscale is used as the background color, and red, green, and blue are used as the overlay colors. Specifically, the LBP peak area is set to red, the Gabor high response area is set to green, and the CD8+ high density area is set to blue, as shown in FIG. Figure 4-6 shown.

[0071] Step S2: Biological Correlation. The 3D fusion image obtained in step S1.2 is calibrated with the multi-center device parameter library and reference standard data, and the correlation between biological markers and texture parameters is obtained through a machine learning model. The details are as follows:

[0072] Step S2.1: Data preprocessing. A multi-center device parameter library is constructed using scanning protocols from different hospitals and devices (e.g., MRI field strength and CT slice thickness). Reference standard data is determined using gold standard pathology results or known biomarkers (e.g., PD-L1 expression and hypoxic area markers). The multi-center device parameter library and reference standard data are then calibrated to obtain a corresponding unified image domain. The details are as follows:

[0073] Step S2.1.1: Parameter library standardization. That is, based on the data values ​​in the multi-center device parameter library and the characteristic mean and standard deviation of the current device / reference device in healthy tissue, the standardized characteristic values ​​corresponding to each data value in the multi-center device parameter library are obtained, specifically:

[0074]

[0075] in: is the standardized eigenvalue, Characteristic values ​​measured for the original equipment, is the characteristic mean of the current device in healthy tissue, is the characteristic standard deviation of the current device in healthy tissue, is the standard deviation of healthy tissue of the reference device, is the mean value of healthy tissue of the reference device.

[0076] In the specific implementation process, the characteristic mean and characteristic standard deviation of the current equipment such as Siemens CT in healthy tissue are 40HU and 25HU respectively, while the healthy tissue mean and healthy tissue standard deviation of the reference equipment such as GE CT are 50HU and 30HU respectively. The original pixel value is 80HU, and the corresponding standardized characteristic value is 98HU.

[0077] Step S2.1.2: Protocol Conversion. This involves constructing a structural model using the liver model (including portal vein grade data) from the 4D XCAD phantom library and the tumor classification ratio provided in the pathology report (e.g., 60% invasive, 30% exophytic, and 10% mixed). It is important to note that during the structural modeling process, a digital voxel matrix containing 5-20 mm tumors is generated for the generated phantom and the respiratory motion curve is embedded, resulting in a customized digital phantom that includes anatomical structure, T1 / T2 parameter maps, elastic modulus distribution, and dynamic respiratory displacement field.

[0078] Furthermore, the sequence parameters of center A (3T) and center B (1.5T) are input into the obtained customized digital phantom, and the corresponding signal intensity distribution matrix is ​​obtained. At the same time, the obtained signal intensity distribution matrix is ​​processed through the constructed loss function to obtain the corresponding calibration slope and calibration intercept, specifically:

[0079]

[0080] in: is the calibration slope, is the original signal strength measured at center B, To calibrate the intercept, is the standard signal reference value of center A, is the L2 norm.

[0081] Furthermore, based on the obtained calibration slope and calibration intercept, and according to the original device signal value, the corresponding standardized signal value is determined, specifically:

[0082]

[0083] in: is the normalized signal value, is the original device signal value, is the calibration slope, is the calibration intercept.

[0084] Step S2.1.3: Image domain unification. The calibration slope and intercept obtained in step S2.1.2, the protocol-converted normalized image, and the dual-center paired DICOM file are used as input to the CycleGAN model. The output is the corresponding synthesized feature image, i.e., the microtexture feature image.

[0085] Step S2.2: Biological Validation. This involves combining the microtexture feature image obtained in step S2.1.3 with the biomarker monitoring data to obtain a fused data matrix. This fused data matrix is ​​then used as input to the random forest model, which outputs the correlation between the corresponding texture parameters and the biomarkers. The details are as follows:

[0086] Step S2.2.1: Spatial transcriptome matching. Using the DAPI-stained image from the Nanostring digital pathology as a fixed reference, align the transcriptome array with the pathology slide and obtain the mean and standard deviation of the image area covered by each transcriptome spatial unit. The spatial variation characteristics of the tumor microenvironment are determined by obtaining the mean and standard deviation. Specifically:

[0087]

[0088] in: is the heterogeneity index, is the characteristic measurement value of the m-th spatial unit, is the characteristic average of all spatial units, is the global standard deviation, is the total number of spatial units analyzed, The spatial unit index.

[0089] In the specific implementation process, there are 5 spatial units, and the characteristic measurement values ​​corresponding to each spatial unit are 1.1*10 -3 mm 2 / s, 0.9*10 -3 mm 2 / s, 1.5*10 -3 mm 2 / s, 0.8*10 -3 mm 2 / s and 1.2*10 -3 mm 2 / s. Therefore, the corresponding characteristic average value is 1.1*10 -3 mm 2 / s, the global standard deviation is 0.25, that is, the corresponding heterogeneity index is -0.2.

[0090] Step S2.2.2: CTC proteomic association analysis. This involves temporally correlating CTC features with imaging markers through Spearman correlation analysis to obtain the corresponding Spearman rank correlation coefficient, specifically:

[0091]

[0092] in: is the Spearman rank correlation coefficient, is the sample size, is the rank difference of the g-th data point.

[0093] During the specific implementation process, a numerical table of CTC counts and ADC values ​​of 10 patients with liver cancer was obtained, as shown in Table 1 below:

[0094] Table 1: CTC count and ADC value table

[0095]

[0096] According to the data in Table 1 above, the corresponding Spearman rank correlation coefficient is 0.31.

[0097] Step S3: Predictive modeling. This involves fusing the microtexture feature image obtained in step S2.1.3 with the real-time monitoring data, and constructing a spatiotemporal heterogeneity-immune clonal diversity joint model. Based on this spatiotemporal heterogeneity-immune clonal diversity joint model, the corresponding risk score is obtained. The details are as follows:

[0098] Step S3.1: Data fusion. This involves fusing the microtexture feature image obtained in step S2.1.3 with the real-time monitoring data. Specifically, in the process of constructing the graph structure, each patient is treated as a central node, and each child node includes an imaging feature node, a liquid biopsy node, and an immune node. At the same time, the connection weights between adjacent nodes are set based on the apparent diffusion coefficient heterogeneity index in the microtexture feature image, specifically:

[0099]

[0100] in: is the connection weight between nodes, is the apparent diffusion coefficient heterogeneity index, is the T cell receptor clonal diversity index, is the clinical critical threshold, is the number of circulating tumor cells, is an operator.

[0101] Furthermore, through the GNN model, the central nodes are connected, and the feature vectors of each node in each layer are determined by setting the connection weights between the nodes. Specifically:

[0102]

[0103] in: is the feature vector of node v in the l+1 layer, is the activation function, is the trainable weight matrix of layer l, For vector splicing, is the feature vector of node v in layer l, is the node index, is the set of neighboring nodes of node v, is the connection weight between nodes u and v.

[0104] Step S3.2: Construct a hierarchical model. That is, based on the treatment response probability predicted by the GNN model in step S3.1, compare it with the actual treatment response to obtain the corresponding target loss function, specifically:

[0105]

[0106] in: is the target loss function, is the true label of the rth sample, is the model-predicted probability of treatment response for the rth sample, is the total number of samples, is the sample index, is the tree structure complexity penalty coefficient, is the total number of leaf nodes, is the L2 regularization coefficient, is the L2 norm squared of the weight.

[0107] During the specific implementation process, there are clinical implementation data of two patients as shown in Table 2 below, as shown in the following table:

[0108] Table 2: Clinical implementation data table

[0109]

[0110] At the same time, the initial predicted treatment response probability corresponding to patient P-011 is 0.6, and the initial predicted treatment response probability corresponding to patient P-012 is 0.4, so the corresponding target loss function is 1.12.

[0111] Step S4: Generate a treatment plan. This involves adjusting the treatment cycle through reinforcement learning based on the target loss function obtained in step S3.2 and real-time surgical navigation data (such as fluorescence imaging and respiratory motion tracking) to obtain a corresponding personalized treatment plan. The details are as follows:

[0112] Step S4.1: Toxicity prediction. Based on the elastic modulus map and ADC map, multiple ROI regions are marked in the tumor core area. Based on the elastic modulus in each ROI region, the corresponding elastic modulus standard deviation is obtained. Simultaneously, based on the obtained elastic modulus standard deviation, the corresponding probability of radiation pneumonitis is determined. Furthermore, based on the determined probability of radiation pneumonitis, toxicity levels are specifically categorized. In this embodiment, there are four toxicity levels, each corresponding to a different probability of radiation pneumonitis. Specifically, when the probability of radiation pneumonitis is 0.4-0.6, the corresponding toxicity level is 1. When the probability of radiation pneumonitis is 0.61-0.8, the corresponding toxicity level is 2. When the probability of radiation pneumonitis is 0.81-0.9, the corresponding toxicity level is 3. When the probability of radiation pneumonitis is greater than 0.9, the corresponding toxicity level is 4.

[0113] In this embodiment, the formula for obtaining the probability of occurrence of radiation pneumonitis is specifically as follows:

[0114]

[0115] in: is the probability of radiation pneumonitis, is the base of natural logarithms, is the standard deviation of the elastic modulus, is the coefficient of variation of the ADC value.

[0116] Furthermore, the obtained probability of radiation pneumonitis is compared with a preset probability threshold range, and the corresponding clinical response is determined based on the comparison result. Specifically:

[0117] If the obtained probability of radiation pneumonitis is less than the lower limit of the preset probability threshold, the current radiotherapy plan is continued. If the obtained probability of radiation pneumonitis is within the preset probability threshold, the frequency of high-resolution CT monitoring is increased. If the obtained probability of radiation pneumonitis is greater than the upper limit of the preset probability threshold, the dose is reduced and serum testing is performed.

[0118] During the specific implementation, the standard deviation of the elastic modulus was obtained to be 2.73, and the coefficient of variation of the ADC value was -0.25, corresponding to a probability of 0.78 for radiation pneumonitis. Furthermore, in this embodiment, the preset probability threshold range is set to (0.4, 0.65). In other words, if the current probability of radiation pneumonitis is greater than the upper limit of the preset probability threshold, the dose should be reduced and serum testing should be performed. The current corresponding toxicity level is 2.

[0119] Step S4.2: Solution optimization. That is, based on the target loss function obtained in step S3.2 and the toxicity level determined in step S4.1, the corresponding reward function is obtained through reinforcement learning strategy, specifically:

[0120]

[0121] in: is the immediate reward value, is the efficacy safety factor, is the target loss function, is the objective response rate weight, is the objective response rate, is the toxicity penalty coefficient, Toxicity level.

[0122] Specifically, based on the instant reward value obtained, the current solution is optimized through the Q-learning algorithm, as follows:

[0123]

[0124] in: is the updated state-action value function, is the state-action value function, is the learning rate, is the immediate reward value, is the discount factor, is the maximum expected value of the next state, is the current state, To perform an action.

[0125] During implementation, the action corresponding to shortening the treatment cycle to 14 days is set to 1, the action corresponding to maintaining the standard 21-day treatment cycle is set to 2, and the action corresponding to extending the treatment cycle to 28 days is set to 3. Specifically, during initial treatment, the action corresponding to all new patients is set to 2. Therefore, when the immediate reward value is greater than 0.7, the state-action value function corresponding to shortening the treatment cycle to 14 days is 0.6, the state-action value function corresponding to extending the treatment cycle to 28 days is 0.5, and the state-action value function corresponding to maintaining the standard 21-day treatment cycle is 0.4. In other words, the treatment plan for this patient is to shorten the treatment cycle to 14 days.

[0126] This embodiment further provides a tumor region image texture analysis system, which uses the above-mentioned tumor region image texture analysis method.

[0127] Example 2

[0128] This embodiment provides a method and system for analyzing texture in tumor region images. The specific implementation method is the same as in Example 1, except that, in step S1.1, the MRI T1-weighted image, elastic modulus map, and CD8+ immunohistochemistry slice image are preprocessed to obtain the corresponding binary mask matrix and coordinate point set. The present invention is described below with reference to the specific implementation methods of this embodiment.

[0129] In this embodiment, the T1-weighted image, elastic modulus map, and CD8+ immunohistochemistry slice image are preprocessed to obtain the corresponding binary mask matrix and coordinate point set, as follows:

[0130] Step S1.1.1: Data normalization. Using the MRI T1-weighted image as the spatial reference, the elastic modulus map and the CD8+ immunohistochemical slice images were registered to the same coordinate system using an affine transformation. The CD8+ immunohistochemical slice images were downsampled to maintain the same voxel size as the MRI T1-weighted image, and the spatial continuity of the cell distribution was preserved using the Anczos interpolation method.

[0131] Furthermore, the 2 mm slice-thick MRI data were linearly interpolated in the Z-axis direction to generate 0.5 mm virtual slices, and the elastic modulus map was resampled based on the original phase data of the MRE T1-weighted image to match the 0.5 mm virtual slices.

[0132] Step S1.1.2: Generate an elastic modulus mask. That is, based on the standardized elastic modulus map obtained in step S1.1.1, compare the modulus values ​​corresponding to each spatial position in the elastic modulus map with a preset modulus threshold, and mark each spatial position in the elastic modulus map based on the comparison result, specifically:

[0133] When the modulus size corresponding to the spatial position is less than the preset modulus threshold, the spatial position is marked as 1. Conversely, when the modulus size corresponding to the spatial position is not less than the preset modulus threshold, the spatial position is marked as 0.

[0134] Specifically, the marked value of each spatial position is the generated elastic modulus value. It is worth noting that the preset modulus threshold in this embodiment is specifically set according to the pathological subtype. Specifically, in the process of setting the preset modulus threshold, at least 200 cases of each subtype confirmed by pathology are classified according to the tumor subtype, and the corresponding elastic modulus distribution is statistically analyzed. At the same time, the 25th percentile value in the statistical elastic modulus distribution is used as the lower threshold of the preset modulus threshold. At the same time, for mixed tumors, the geometric mean is used as the corresponding preset modulus threshold.

[0135] In the specific implementation, the elastic modulus of hepatocellular carcinoma is 5 kPa, the elastic modulus of cholangiocarcinoma is 6.5 kPa, and the elastic modulus of metastatic tumors is 7 kPa. Therefore, the preset modulus threshold in this embodiment is set to 5 kPa. In other words, when the modulus corresponding to a spatial position in the elastic modulus map is less than 5 kPa, the spatial position is marked as 1, and when the modulus corresponding to a spatial position in the elastic modulus map is not less than 5 kPa, the spatial position is marked as 0.

[0136] Step S1.1.3: CD8+ immunolabeling. The CD8+ cell localization coordinates in the CD8+ immunohistochemistry slice image are converted to a global coordinate system aligned with the MRI T1-weighted image, and the corresponding nuclear density is obtained, specifically:

[0137]

[0138] in: For coordinates The cell density at is the bandwidth parameter, is the target location where the density is to be estimated, is a non-negative function for weight distribution, is the total number of cell events in the current analysis window, is the center coordinate of the ith immune cell, An index of immune cells.

[0139] Furthermore, the obtained nuclear density is compared with the preset cell density threshold, and the center coordinates of the immune cells are deleted based on the comparison result. Specifically:

[0140] When the obtained nuclear density is less than the preset cell density threshold, the center coordinates of the immune cell corresponding to that nuclear density are deleted. Conversely, when the obtained nuclear density is not less than the preset cell density threshold, the center coordinates of the immune cell corresponding to that nuclear density are retained. Specifically, in this embodiment, the immune cell corresponding to the retained center coordinates is a CD8+ immune-labeled cell.

[0141] Although embodiments of the present invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and variations may be made to these embodiments without departing from the principles and spirit of the invention, and that the scope of the invention is limited by the accompanying embodiments and their equivalents.

Claims

1. A method for texture analysis of tumor region images, characterized in that: Includes: S1: Image preprocessing: preprocessing the acquired multimodal data to obtain preprocessed multimodal data, and performing feature extraction on the preprocessed multimodal data to obtain a feature image; S2: Biological Correlation: Calibrate the feature image with a multi-center device parameter library and reference standard data, and obtain the correlation between biological markers and texture parameters through a machine learning model, including: S2.1: Data preprocessing: Calibrate the constructed multi-center equipment parameter library and reference standard data to obtain a unified image domain; S2.2: Biological Validation: Combining the unified image domain with the biomarker monitoring data to obtain a fused data matrix, and using the fused data matrix as input to a random forest model to output correlations between texture parameters and biomarkers; S3: Predictive modeling: Based on the unified image domain and real-time monitoring data, a spatiotemporal heterogeneity-immune clonal diversity joint model is constructed to obtain a risk score; S4: Generate treatment plan: Based on the risk score and real-time surgical navigation data, adjust the treatment cycle through reinforcement learning to obtain a personalized treatment plan.

2. The method for analyzing tumor region image texture according to claim 1, characterized in that: The feature image is obtained, including: S1.1: Image processing: Align the T1-weighted image, elastic modulus map, and CD8+ immunohistochemical slice image to the same coordinate system using predefined anatomical landmarks. Resample the voxel sizes of the elastic modulus map and CD8+ immunohistochemical slice image to match those of the T1-weighted image. Obtain a binary mask matrix for the elastic modulus map and an immune heat map for the CD8+ immunohistochemical slice image. S1.2: Microtexture extraction: Based on the T1-weighted image, binary mask matrix and immune heat map, a three-dimensional multimodal data cube is constructed, and the three-dimensional multimodal data cube is processed using a gray-level co-occurrence matrix, local binary pattern and Gabor filter to obtain a multimodal feature matrix and a three-dimensional fusion image.

3. The method for analyzing tumor region image texture according to claim 2, characterized in that: The three-dimensional fusion image includes a background color and an overlay color. The background color is set to the grayscale of the T1 weighted image, and the overlay color is set to different colors according to the LBP peak area, the Gabor high response area and the CD8+ high density area.

4. The method for analyzing tumor region image texture according to claim 2, characterized in that: Obtaining the binary mask matrix of the elastic modulus map and the immune heat map of the CD8+ immunohistochemistry slice image includes: S1.1.1: Data normalization: Based on the T1-weighted images, register the elastic modulus map and CD8+ immunohistochemistry slice images to the same coordinate system using affine transformation, and downsample the CD8+ immunohistochemistry slice images; S1.1.2: Generate elastic modulus mask: Compare the magnitude of each modulus in the standardized elastic modulus map with a preset modulus threshold, and mark each spatial position in the standardized elastic modulus map based on the comparison result, specifically: When the modulus at the spatial position is smaller than the preset modulus threshold, the spatial position is marked as 1; otherwise, the spatial position is marked as 0; S1.1.3: CD8+ immune markers: Based on the CD8+ cell localization coordinates and the global coordinate system of the T1-weighted image, obtain the nuclear density size, compare the nuclear density size with the preset cell density threshold, and delete the center coordinates of the immune cells based on the comparison results, specifically: When the nuclear density is smaller than the preset cell density threshold, the central coordinates of the immune cells corresponding to the nuclear density are deleted; otherwise, the central coordinates of the immune cells corresponding to the nuclear density are retained.

5. The method for analyzing tumor region image texture according to claim 1, characterized in that: The unified image domain is obtained, including: S2.1.1: Parameter library standardization: Obtain standardized characteristic values ​​of the multi-center device parameter library based on the data values ​​in the multi-center device parameter library and the characteristic means and standard deviations of the current device / reference device in healthy tissue; S2.1.2: Protocol Conversion: Perform structural modeling based on the liver module and tumor classification ratio to obtain a customized digital phantom. Simultaneously, obtain the calibration slope and intercept using the sequence parameters of Center A and Center B. Based on the calibration slope and intercept, determine the standardized signal value, specifically: in: is the normalized signal value, is the original device signal value, is the calibration slope, is the calibration intercept; S2.1.3: Image domain unification: The standardized signal values ​​and the dual-center paired DICOM files are used as inputs to the CycleGAN model, and the synthesized feature images are obtained as output.

6. The method for analyzing tumor region image texture according to claim 1, characterized in that: The outputs obtained were correlations between texture parameters and biological markers, including: S2.2.1: Spatial transcriptome matching: Based on the DAPI staining image, the transcriptome dot matrix is ​​aligned with the pathological section, and the mean and standard deviation of the image area covered by the transcriptome spatial unit are obtained to determine the spatial variation characteristics of the tumor microenvironment. in: is the heterogeneity index, is the characteristic measurement value of the m-th spatial unit, is the characteristic average of all spatial units, is the global standard deviation, is the total number of spatial units analyzed, is the spatial unit index; S2.2.2: CTC proteomic association analysis: Spearman correlation analysis was used to temporally correlate CTC characteristics with imaging markers to obtain the Spearman rank correlation coefficient.

7. The method for analyzing tumor region image texture according to claim 1, characterized in that: Obtain risk scores, including: S3.1: Data Fusion: Each patient is set as a central node, and the connection weights between adjacent nodes are set based on the apparent diffusion coefficient heterogeneity index of the unified image domain. At the same time, the set central nodes are connected through the GNN model, and the feature vectors of each node at each layer are determined based on the connection weights; S3.2: Constructing a hierarchical model: Based on the treatment response probability of the GNN model and the actual treatment response, the target loss function is obtained, specifically: in: is the target loss function, is the true label of the rth sample, is the model-predicted probability of treatment response for the rth sample, is the total number of samples, is the sample index, is the tree structure complexity penalty coefficient, is the total number of leaf nodes, is the L2 regularization coefficient, is the L2 norm squared of the weight.

8. The method for analyzing tumor region image texture according to claim 1, characterized in that: Get a personalized treatment plan, including: S4.1: Toxicity prediction: Based on the elastic modulus graph and ADC graph, obtain the elastic modulus standard deviation, and determine the probability of radiation pneumonitis based on the elastic modulus standard deviation. At the same time, classify the toxicity level based on the probability of radiation pneumonitis; S4.2: Solution optimization: Based on the risk score and toxicity level, a reward function is obtained through a reinforcement learning strategy, and the current solution is optimized based on the reward function and the Q-learning algorithm.

9. The method for analyzing tumor region image texture according to claim 8, characterized in that: The probability of occurrence of radiation pneumonitis is compared with a preset probability threshold range, and the corresponding clinical response is determined based on the comparison result, specifically: When the probability of radiation pneumonitis is less than the lower limit of the preset probability threshold, the current radiotherapy plan is continued; when the probability of radiation pneumonitis is within the preset probability threshold, the frequency of high-resolution CT monitoring is increased; when the probability of radiation pneumonitis is greater than the upper limit of the preset probability threshold, the dose is reduced and serum testing is performed.

10. A tumor region image texture analysis system, characterized in that: The tumor region image texture analysis method according to any one of claims 1 to 9 is used.

Images