Microenvironment formed by pig muscle fibers and interactive characteristic quantification method thereof

Abstract

The invention discloses a microenvironment formed by pig muscle fibers and an interactive characteristic quantification method thereof. The method comprises the following steps: S1, extracting image information; s2, diffusion circle microenvironment extraction; s3, connected domain microenvironment extraction; s4, muscle fiber interaction characteristic measurement; s5, microenvironment characteristic index measurement; s6, microenvironment classification: similarity clustering is carried out on the diffusion circle and the connected domain according to the rotation overlapping degree, the connected domain microenvironment feature vector representation, weighted k-means and PCAk-means clustering; and S7, measuring the change trend of the diffusion circle characteristic along with the increase of the radius. Compared with the prior art, the method can provide quantification of apparent characteristics, which are not considered yet, of pig muscle fibers, the characteristics of the microenvironment where the muscle fibers are located, the interaction characteristics between the muscle fibers and the position characteristics of the muscle fibers, and the method is set up to recognize typical microenvironment categories based on images.

Description

Technical Field

[0001] This application belongs to the field of cell image information extraction technology, specifically involving a method for quantifying the microenvironment and interactive features of porcine muscle fibers. Background Technology

[0002] As one of my country's main meat sources, skeletal muscle is the primary source of meat production, and muscle fibers are a crucial component. Extracting the phenotypic shape of microscopic muscle fibers and further investigating the correlation between these shapes and macroscopic indicators of skeletal muscle growth and development, as well as meat quality, helps to elucidate the mechanisms of meat formation and the phenotypic changes in muscle fibers during growth, development, wound repair, and environmental interactions. Pigs also play an important role as model organisms: pig fat types can be directly correlated with human fat types, and recent significant progress in the transplantation of pig hearts and kidneys into humans demonstrates that biological analysis of pigs contributes to understanding human biological mechanisms, and that modifying pigs can create valuable medical components for humans. Research on pig skeletal muscle growth, development, wound repair, and environmental interactions at the muscle fiber level contributes to the study of corresponding human muscle components. In conclusion, research on pig muscle fiber phenotypes has significant relevance, with the extensive extraction of phenotypes being fundamental.

[0003] Cells do not function independently, but work together with other cells to regulate the physiological activities of animals. A unit with a specific function composed of a group of cells in a specific way is called a microenvironment or a repetitive unit. The characteristic of the microenvironment is usually the repetition of local areas with similar cell composition, and their identification can be carried out based on this characteristic. The identification of the microenvironment helps to study and understand the synergistic mechanism of cells and the hierarchical understanding of functional units in biological tissues. Since it is known that the epigenetic characteristics and gene expression characteristics of cells are affected by the regulation of other cells (called cell communication / regulation), the degree and characteristics of this regulation are obviously affected by the spatial distribution, contact mode and contact degree of cells, which is called the spatial interaction characteristics of cells. For example, the literature [1] used 10×Genomics spatial RNA sequencing technology, muscle tissue imaging method combined with pseudo-temporal analysis and other techniques to analyze how human muscle tissue granulomas affect the surrounding muscle tissue at short and long distances, revealing the gradient of inflammation and fibrosis near the granuloma. Reference [2] combined a large amount of single-cell transcriptome data with high-throughput spatial transcriptome data to generate a mouse muscle repair trajectory, discovering key events and spatiotemporal characteristics of muscle repair at the molecular level, and identifying cell-cell dynamic interactions. Reference [3] used 10×Genomics spatial RNA sequencing technology and single-cell transcriptomics technology to generate a high-resolution molecular atlas of muscle tissue in malnourished mice, clearly defining the cell types and gene spatiotemporal expression patterns associated with different regeneration stages in the mouse model, discovering important spatiotemporal features of inflammation and fibrosis at the cellular and molecular levels, and generating a spatial database of the muscle tissue microenvironment to explore how local cell-cell communication promotes the progression of related diseases, such as exploring specific cell types or activated signaling pathways in damaged or regenerated areas, and providing information on intercellular crosstalk-mediated diseases. Quantifying the interaction characteristics of muscle fibers enables correlation analysis with other properties of muscle fibers (such as gene expression information obtained from spatial transcriptomics and cellular biochemical properties, as mentioned in the literature) in skeletal muscle research. This allows for the elucidation of spatial regulatory mechanisms among muscle fibers, which is of great significance for understanding skeletal muscle growth, repair, and fattening. It can also help improve the efficiency of pig farming and provide molecular-level inspiration for human skeletal muscle research from the perspective of model animals. This invention proposes a quantitative scheme for the phenotypic characteristics of porcine muscle fibers, focusing primarily on the quantification of the microenvironment and interaction characteristics of muscle fiber pairs.

[0004] Early research on muscle fiber phenotypes primarily focused on the proportion of different muscle fiber types and conventional size indicators, such as area and diameter, analyzing the trends of these indicators with physiological location and growth and development. In the 1980s and 1990s, research on the geometric distribution patterns of muscle fibers mainly involved fitting planar images of muscle fibers with regular polygonal grids to construct fiber adjacency graphs (FAGs), then defining specific clusters and establishing statistics to determine the degree of geometric, random, and clustered distributions. Different studies used different grid fitting methods, such as hexagonal grid fitting, Gabriel plots and Delaunay triangulation, radial axis triangulation and generalized Dirichlet mosaic, as well as different statistics, such as average cluster size, number of closed fibers, and co-dispersion index. Additionally, some statistics based on muscle fiber distance were also used. These studies corrected the perception that "muscle fiber types in healthy muscle tissue tend to be randomly distributed." Subsequent research has gradually begun to study the distribution patterns of muscle fibers from the perspective of function fitting, including the use of binary Markov random fields and the recently proposed generalized additive function fitting. These existing methods and techniques have the following characteristics: (1) they do not deeply address the identification and feature quantification of the muscle fiber microenvironment; (2) they do not deeply address the quantification of interaction features between muscle fibers; (3) most methods are based on planar segmentation data rather than precise muscle fiber contours; (4) existing methods do not consider the heterogeneous nucleus situation of muscle fibers, and classify muscle fibers as absolute membership rather than membership degrees of 0 to 1. This invention will focus on the above problems and propose a scheme for the identification and quantification of the muscle fiber microenvironment and muscle fiber interaction features based on precise muscle fiber contours.

[0005] This approach establishes a system based on the type and contour data of muscle fibers in stained porcine muscle tissue images to mine local circular microenvironments with similar cellular components and geometric distributions, as well as connected regions of similar cells with similar characteristics. These spatial microenvironments help identify specific functional modules and cellular synergies. Based on this, the invention designs an algorithm system to measure the interaction characteristics between muscle fibers according to their type and contour data. This algorithm can quantify the degree and characteristics of spatial interaction between each muscle fiber and its neighboring or contiguous microenvironments, as well as its relative centrality in various local environments. This method can also be used for feature analysis of other non-overlapping cell contours.

[0006] [1]H.Lequain,C.Dégletagne,N.Streichenberger,J.Valantin,T.Simonet,L.Schaeffer,P.Sève,P.Leblanc.Spatial Transcriptomics Reveals Signatures ofHistopathological Changes in Muscular Sarcoidosis[J].Cells,2023,12(23):2747.

[0007] [2]DWMckellar,LDWalter,LTSong,M.Mantri,MFZWang,I.DeVlaminck,BDCosgrove.Large-scale integration of single-cell transcriptomicdata captures transitional progenitor states in mouse skeletal muscleregeneration[J].Communications Biology,2021,4(1):1280.

[0008] [3]MJStec,Q.Su,C.Adler,L.Zhang,DRGolann,NPKhan,P.Lampros,S.Panagis,A.Villalta,M.Ni,Y.Wei,JRWalls,AJMurphy,GDYancopoulos,GSAtwal,S.Kleiner,G.Halasz,MWSleeman.A cellular and molecular spatial atlas of dystrophic muscle[J].Proceedings of the National Academy ofSciences,2023,120(29):e2221249120. Summary of the Invention

[0009] This application provides a method for quantifying the microenvironment and interaction characteristics of porcine muscle fibers, aiming to solve the problems of how to effectively identify the spatial microenvironment of muscle fibers in porcine muscle tissue and measure its characteristics, as well as the interaction characteristics between muscle fibers.

[0010] A method for quantifying the microenvironment and interaction characteristics of porcine muscle fibers, the method comprising the following steps:

[0011] S1: Image information extraction. Based on the cross-sectional imaging of muscle fibers with pre-existing muscle fiber ID, category, contour coordinate set, and geometric center coordinates, calculate the center distance matrix and shortest distance matrix between muscle fibers, and use the minimum bounding rectangle to calculate the short axis data of each muscle fiber to construct the adjacency matrix.

[0012] S2: Diffusion circle microenvironment extraction: Starting from the geometric center of each muscle fiber, three diffusion circles are extracted with radii of 1.5, 2.5, and 3.5 times the average diameter of all muscle fibers, respectively. Starting from the geometric center of each muscle fiber, the adjacency relationship of muscle fibers is identified based on the distance matrix of the muscle fibers. Three adjacent diffusion circles are generated with the first adjacent layer, the first and second adjacent layers, and the first, second, and third adjacent layers as the range.

[0013] S3: Connected domain microenvironment extraction, based on the distance matrix of muscle fibers to identify connected domains formed by the adjacency relationship of the same type of muscle fibers;

[0014] S4: Muscle fiber interaction feature measurement, which measures the degree of contact between a muscle fiber and its various types of muscle fibers, measures its multi-layer adjacency information in the diffusion circle, and its centrality in the connected domain.

[0015] S5: Microenvironment characteristic indicators, including indicators of the size and shape of connected domains, indicators of the composition ratio and distribution pattern of various types of muscle fibers in the diffusion circle, and statistical measures of muscle fiber properties in these two microenvironments.

[0016] S6: Microenvironment classification. Based on the degree of rotational overlap and the weighted Euclidean distance between microenvironment feature vectors, weighted k-means and PCAk-means clustering are used to perform similarity clustering on the diffusion circle and connected components.

[0017] S7: Measurement of the changing trend of diffusion circle characteristics with increasing radius. For each muscle fiber, a diffusion circle is formed at 1.5, 2.5, and 3.5 times the average diameter of all muscle fibers, and the first adjacent layer, the first and second adjacent layers, and the first, second, and third adjacent layers are used as the range. The proportional vector index in the three adjacent diffusion circles is generated and the first and second differences are calculated.

[0018] Optionally, the distance matrix construction in S1 includes the following operations:

[0019] Calculate the center distance and shortest distance between any two muscle fibers in the image, where the center distance between the i-th muscle fiber and the j-th muscle fiber is denoted by cd. ij express:

[0020]

[0021] Shortest distance is MD ij express:

[0022]

[0023] Where Ω i and Ω j Let CD represent the sets of contour coordinates of the i-th and j-th muscle fibers, respectively, and let CD = cd. ij Distance matrix MD = {md ij}

[0024] Optionally, the extraction of the diffusion circle microenvironment in S2 includes the following operations:

[0025] S2.1: Construct an adjacency matrix for all muscle fibers. Initialize the adjacency matrix ADJ to the distance matrix MD. Denote the elements whose values ​​are less than 1 / 2 of the minimum length of the minor axis of all muscle fibers as MSD and set their values ​​to 1. This completes the construction of the adjacency matrix ADJ.

[0026] The adjacency matrix ADJ has elements that satisfy:

[0027]

[0028] Where md ij The shortest distance matrix for muscle fibers is MD = {md} ij The elements of} represent the shortest distance between the i-th and j-th muscle fiber contours, and MSD is the minimum minor axis of all muscle fibers in the image;

[0029] S2.2: The first type of diffusion circle is uniquely represented by the center point of the microenvironment and the absolute radius of the microenvironment. The center point of the microenvironment is the geometric center of a certain muscle fiber. Match the muscle fiber ID, and at the same time, exclude circles whose contours exceed the image boundary.

[0030] S2.3: The second type of diffusion circle is uniquely represented by the center point of the microenvironment and the adjacency radius of the microenvironment. The adjacency radius is defined by the distance between the muscle fibers in the adjacency matrix. If muscle fiber i needs to move k steps to reach muscle fiber j through adjacency relationship, then their distance is k. For each valid diffusion circle with three radius lengths of adjacency distances of 1, 2, and 3, the muscle fiber ID where the center is located, the coordinates of the center, and the adjacency radius of the diffusion circle are obtained as the identifier of this microenvironment. The muscle fiber IDs contained therein are calculated. Taking the diffusion circle of the i-th muscle fiber with an adjacency radius of k as an example, the calculation method of the muscle fiber IDs contained therein is as follows:

[0031] find(ADJ(i,:)+ADJ 2(i,:)+…+ADJ k (i,:))

[0032] ADJ k (i,:) represents the i-th row of the adjacency matrix raised to the power of k, and the find function returns the position of the non-zero element in it.

[0033] Optionally, the extraction of the connected domain microenvironment in S3 includes the following operations:

[0034] S3.1: Construct the adjacency matrix ADJ for all muscle fibers of class k. (k) ADJ (k) =ADJ(N k N k ), where N k Let be the ID of all k-th type muscle fibers in the image;

[0035] S3.2: Connected Component Identification: For each type of muscle fiber, a breadth-first search algorithm is used to identify connected regions.

[0036] Optionally, the measurement of muscle fiber interaction characteristics in S4 includes the following operations:

[0037] S4.1: Calculate the ratio of the contact circumference of each muscle fiber to that of each type of muscle fiber;

[0038] S4.2: Calculate the weighted sum of the contact ratios between each muscle fiber and surrounding muscle fibers, and select specific weights to reflect the interaction strength at a specific level;

[0039] S4.3: Calculate the information within the circle radiating outward from the muscle fiber as the center. Taking the k-th muscle fiber as an example, the information of the radiating circle radiating outward from the k-th muscle fiber as the center is condensed into a matrix consisting of distance, angle, category, and area occupied.

[0040] S4.4: Measures the interaction strength between muscle fibers and various types of muscle fibers in a specified local area;

[0041] S4.5: Measures the degree of centrality of a muscle fiber within a connected domain formed by similar muscle fibers.

[0042] Optionally, the S5 microenvironment characteristic index measurement includes the following operations:

[0043] S5.1: Quantification of diffusion circle microenvironment characteristics, including:

[0044] S5.1.1: The proportion of the total area of ​​each type of muscle fiber within the diffusion circle. If the obtained value is a category membership degree, it is classified according to a threshold of 0.5. Taking the diffusion circle of the i-th muscle fiber with radius or adjacent concentric circles of r as an example, it is denoted as ARatio. i r ;

[0045] S5.1.2: The proportion of arc occupied by each type of muscle fiber within the diffusion circle;

[0046] S5.1.3: The adjacency information matrix of the i-th muscle fiber radius or adjacency layer r constructed above.

[0047] S5.1.4: Statistical values ​​of interaction attributes of various types of muscle fibers;

[0048] S5.2: Calculate the autocorrelation Moran coefficient of the muscle fiber index values ​​within the diffusion circle;

[0049] S5.3: Calculate the local distance matrix entropy and included angle matrix entropy values ​​of various types of muscle fibers within the diffusion circle;

[0050] S5.4: Measures the degree of aggregation of different types of muscle fibers within the diffusion circle;

[0051] S5.5: Quantification of microenvironment features of connected domains, including the following steps:

[0052] S5.5.1: Quantitative size indicators: area and number of muscle fibers;

[0053] S5.5.2: Quantitative Shape Indicators: Aspect Ratio, Roundness, Fill Level, Rectangularity

[0054] S5.5.3: Calculate the cluster coefficient.

[0055] Optionally, the S6 microenvironment classification includes the following operations:

[0056] S6.1: Use sparse matrix clustering for absolute diffusion circles;

[0057] S6.2: Use weighted k-means and PCA-Kmeans clustering for connected components;

[0058] S6.3: For valid adjacent diffusion circles, the above PCA-Kmeans clustering is used. The selected features are the proportion of the total area of ​​each type of muscle fiber, the proportion of the arc occupied by each type of cell, the entropy values ​​of the local distance matrix and the included angle matrix of each type of muscle fiber, and the degree of aggregation of each type of muscle fiber. These are stretched into a feature vector instead of a matrix. If the adjacency radius is greater than 1, for example, k, then the adjacent diffusion circle features are calculated for the adjacency radii i = 1, 2, ..., k respectively. These features are then concatenated into a feature vector. Next, based on these feature vectors, the above PCA-Kmeans clustering is performed.

[0059] Optionally, the measurement of the change trend of the diffusion circle feature with increasing radius or adjacent concentric layers in S7 includes the following operations:

[0060] S7.1: For a certain muscle fiber i that can generate an effective diffusion circle, for the geometric center of the muscle fiber (x i ,y i Generate diffusion circles with three radii: r1 = 1.5·D avg r² = 2.5·D avg r3 = 3.5·D avg , where D avg This represents the global average equivalent diameter of the current type of muscle fiber;

[0061] S7.2: Using the previously calculated proportional vectors, combine them into a proportional matrix:

[0062]

[0063] S7.3: First-order one-step difference calculation: Taking the k-th type of muscle fiber as an example, using three points... To calculate the second-order difference, first calculate the one-step first-order difference between adjacent points:

[0064]

[0065] S7.4: The second difference is the difference of the first difference:

[0066] S7.5: First-order two-step difference calculation:

[0067] Optionally, the weighted k-means includes the following steps:

[0068] (1) Data Reading and Feature Normalization: The pandas library is used to read the data from the document storing connected component information. For the features that need to be clustered, features = ['circularity', 'aspect ratio', 'fill degree', 'rectangularity', 'total area', 'cluster coefficient'], MinMaxScaler is used for normalization. The normalization formula is:

[0069]

[0070] Where x is the original feature value, that is, the specific value of a certain sample in the dataset for a certain feature;

[0071] (2) Weight settings: Custom weight list w * =[0.1, 0.1, 0.1, 0.1, 0.5, 0.1], because the total area is more intuitive, so the weight is set to the maximum;

[0072] (3) The elbow rule is used to determine the number of clusters, optimal_clusters, and the maximum number of clusters, num_clusters, is set to 20;

[0073] For each k cluster Solve for the objective function:

[0074]

[0075] Where, c1,…, Representing different clustering partitioning schemes, the goal is to find the optimal partitioning scheme that minimizes the objective function value. Double summation: the outermost layer... Summing over all clusters, intermediate layer Sum the squared distances between all samples within each cluster and the cluster center, then calculate the weighted Euclidean distance:

[0076]

[0077] in, Let u be the normalized eigenvector of the j-th connected component. i Let be the center (mean) of the eigenvectors of the connected components in the i-th cluster. The weight of the k-th feature;

[0078] Calculate the second derivative:

[0079]

[0080] SSE″(k)≈SSE(k+1)-2SSE(k)+SSE(k-1),

[0081] Where k∈[2,num_clusters-1], k-1=0 cannot be calculated when k=1; k+1=num_clusters cannot be calculated when k=num_clusters (it exceeds the initially set range);

[0082] Optimal number of clusters:

[0083]

[0084] Optionally, the PCA-Kmeans clustering includes the following steps:

[0085] (1) Data Reading and Feature Selection: Taking connected components as an example, read the connected component data, select 6 original features for analysis: features = ['circularity', 'aspect ratio', 'fill degree', 'rectangularity', 'total area', 'clustering coefficient'], construct the original feature matrix X, let the number of connected components of an image be n, the number of features be p = 6, and the original feature matrix X = [x ij ]∈R n×p , where x ij This represents the j-th eigenvalue of the i-th connected component;

[0086] (2) Z-score standardization: Standardize the original features;

[0087] Among them

[0088] (3) Principal component analysis:

[0089] First, calculate the eigenvalues and eigenvectors of the covariance matrix of the standardized data;

[0090] Covariance matrix:

[0091] Calculation of eigenvalues and eigenvectors: Solve ∑v k = λ k v k , where λ k is the eigenvalue, and v k is the corresponding unit eigenvector;

[0092] Then, perform the calculation of the explained variance ratio. Proportion of the explained variance by a single principal component:

[0093]

[0094] Among them, λ k is the eigenvalue corresponding to the k-th principal component (the k-th eigenvalue), p is the total number of original features, Sum of all principal component eigenvalues;

[0095] Cumulative explained variance ratio:

[0096] where k represents the first k principal components (k ≤ p);

[0097] The cumulative explained variance ratio represents the proportion of the total variance of the original data jointly explained by the first k principal components;

[0098] Determine the number of principal components. When n ≥ p, select k such that the cumulative variance ratio ≥ 85%; when n < p, select k = min(n - 1, p) and the cumulative variance ratio ≥ 80%;

[0099] Calculate the principal component scores. Project the standardized data into the principal component space to obtain the principal component score matrix Z after dimensionality reduction: Z = X scaled ·V, where V = [v 1, v 2, [[ID=6�]]v 3, …, v k ∈ R 0000042;

[0100] (4) Directly use the principal component scores in Z for clustering:

[0101] The elbow rule determines the number of clusters, optimal_clusters, and the maximum number of clusters, num_clusters, is set to 20;

[0102] For each k cluster Solve for the objective function:

[0103]

[0104] To determine the clustering scheme, where C i Let u be the set of indices contained in the i-th cluster. i It is the center vector of the Z matrix, which is the mean of the corresponding principal component score vector;

[0105] Calculate the second derivative:

[0106]

[0107] SSE″(k)≈SSE(k+1)-2SSE(k)+SSE(k-1),

[0108] Where k∈[2,num_clusters-1], k-1=0 cannot be calculated when k=1; k+1=num_clusters cannot be calculated when k=num_clusters (it exceeds the initially set range);

[0109] Optimal number of clusters:

[0110]

[0111] Compared with the prior art, this application has at least the following beneficial effects:

[0112] Compared with existing technologies, this invention can provide quantification of previously unconsidered phenotypic features (biologically significant) of porcine muscle fibers, mainly the microenvironmental features of the muscle fibers, the interaction features between muscle fibers, and the positional features of the muscle fibers. The method is based on image recognition of typical microenvironment categories.

[0113] Compared to existing technologies, the micro-features of the diffusion circle microenvironment provided by this invention include the area ratio of each type of cell, the radian ratio they occupy, the entropy values ​​of their mutual distribution distance and included angle, the average distance ratio of the same type, and the first and second differences of these ratios as the area changes. These features characterize the distribution of other cells around the cells and the trend of this distribution changing with distance. Among them, the Moran index, calculated within the diffusion circle using indicators such as cell type membership and area, characterizes the spatial autocorrelation of these indicators within the microenvironment.

[0114] The microenvironmental characteristics of connected domains of the same type of muscle fibers provided by this invention include aspect ratio, roundness, rectangularity, filling degree, and clustering coefficient, which can quantify the shape pattern and aggregation degree of adjacent clusters of various types of muscle fibers.

[0115] Compared to existing technologies, the out-in contact ratio between muscle fibers measured in this invention can quantify the degree of contact between each muscle fiber and various types of muscle fibers, which may be related to the intercellular communication bandwidth; the measured spatial interaction strength between muscle fibers and various types of muscle fibers characterizes the degree of influence of other muscle fibers at the spatial level; the measured weighted centrality of muscle fibers in connected domains measures their positional characteristics in similar connected regions, which may reflect their niche in similar regulatory chains.

[0116] Compared to existing technologies, the microenvironment and cell interaction characteristics provided in this approach can be correlated with other properties of muscle fibers (such as category, gene expression, and biochemical properties) in skeletal muscle research to elucidate the spatial regulatory mechanisms between muscle fibers. The provided microenvironment clustering methods and results help identify specific functional modules and their constituent structures, aiding in the analysis of functional formation mechanisms. This is of great significance for understanding skeletal muscle growth, repair, and hypertrophy, and can also provide inspiration for human skeletal muscle research from the perspective of model organisms.

[0117] The method developed in this invention has been tested for phenotypic analysis of skeletal muscle fiber distribution in Rongchang pigs and Large White pigs, leading to new discoveries on muscle fiber distribution characteristics between different parts and breeds. Attached Figure Description

[0118] Figure 1 This is an example of the difference in the proportion of arc occupied by type I muscle fibers in the diffusion circle microenvironment as defined in this invention;

[0119] Figure 2 This is a visualization of the connected component feature indicators and diffusion circles defined in this invention.

[0120] Figure 3 This is a framework diagram of the method for identifying muscle fiber microenvironments, measuring microenvironment characteristics, and measuring muscle fiber interaction characteristics according to the present invention.

[0121] Figure 4 This is a flowchart illustrating the process and method framework for measuring and clustering the microenvironmental characteristics of connected domains of similar muscle fibers in this invention.

[0122] Figure 5 This is a flowchart illustrating the process and method framework for measuring and clustering the microenvironmental characteristics of each diffusion circle in this invention.

[0123] Figure 6 This is a flowchart illustrating the measurement process and method framework for muscle fiber interaction indices in this invention.

[0124] Figure 7 This is a clustering test case of the muscle fiber diffusion loop and similar muscle fiber connected domains designed in this invention;

[0125] Figure 8 The method designed in this invention was applied to the differences in muscle fiber indices between Large White pigs and Rongchang pigs.

[0126] Figure 9 The method designed for this invention was applied to the differences in muscle fiber connectivity indices between body parts of Large White pigs and Rongchang pigs.

[0127] Figure 10 The method designed for this invention was applied to the differences in muscle fiber interaction characteristics between body parts of Large White pigs and Rongchang pigs.

[0128] Figure 11 The method designed for this invention is applied to the differences in skeletal muscle regions and between skeletal muscle regions found in Large White pigs and Rongchang pigs;

[0129] Figure 12 The method designed in this invention was applied to the differences in connectivity indicators between the same and different skeletal muscle regions in Large White pigs and Rongchang pigs.

[0130] Figure 13 The method designed in this invention was applied to the differences in muscle fiber interaction indices between the same and different skeletal muscle regions in Large White pigs and Rongchang pigs.

[0131] Figure 14 The method designed for this invention was applied to the Pearson correlation coefficient of muscle fiber area and interaction index in Large White pigs and Rongchang pigs.

[0132] Figure 15 The Spearman correlation coefficient between the type I muscle fiber connectivity feature index and meat quality index was found when the method designed in this invention was applied to the longissimus dorsi muscle of Large White pigs.

[0133] Figure 16 The Spearman correlation coefficient between the spatial interaction characteristics of type I muscle fibers and meat quality indicators found by applying the method designed in this invention to the longissimus dorsi muscle of Large White pigs;

[0134] Figure 17 The Spearman correlation coefficient between the type I muscle fiber connectivity feature index and intramuscular fat content index was obtained by applying the method designed for this invention to the longissimus dorsi and psoas major muscles of Large White and Landrace pigs.

[0135] Figure 18 The Spearman correlation coefficient of the spatial interaction characteristics of type I muscle fibers and intramuscular fat content in the longissimus dorsi and psoas major muscles of Large White and Landrace pigs was obtained by applying the method designed in this invention to them.

[0136] Figure 19 This invention defines and structures the microenvironment characteristic indicators and muscle fiber interaction characteristic indicators.

[0137] Figure 20 The results are for the connected component feature index.

[0138] Figure 21 The results are calculated for the connectivity centrality and diffusion circle partial characteristic indices.

[0139] Figure 22 This is a diagram showing the English abbreviations for skeletal muscles.

[0140] Figure 23 This is a diagram showing the abbreviations of muscle fiber characteristic indicators. Detailed Implementation

[0141] To make the objectives, technical solutions, and advantages of this application clearer, the following detailed description is provided in conjunction with the accompanying drawings and embodiments.

[0142] The method for quantifying the microenvironment and its interaction characteristics of porcine muscle fibers provided in this application includes the following steps:

[0143] S1: Image Information Extraction: Based on the existing cross-sectional imaging of muscle fibers with [muscle fiber ID, category, contour coordinate set, geometric center coordinates], the center distance matrix and the shortest distance matrix between muscle fibers are calculated, and the short axis data of each muscle fiber is calculated using the minimum bounding rectangle to construct the adjacency matrix. Note: (1) How to prepare and label the cross-sectional images of muscle fibers is existing technology and will not be described in detail in this solution; (2) Category labeling can be a definite category attribute or a category membership degree (e.g., the membership degree between type II and type I of muscle fibers obtained based on information such as staining depth); (3) The geometric center can be omitted and can be easily calculated based on the contour data.

[0144] S2: Diffusion Circle Microenvironment Extraction: Starting from the geometric center of each muscle fiber, three diffusion circles are extracted with radii of 1.5, 2.5, and 3.5 times the average diameter of all muscle fibers, respectively. Starting from the geometric center of each muscle fiber, the adjacency relationship of the muscle fibers is identified based on the distance matrix of the muscle fibers. Three adjacent diffusion circles are generated with the first adjacent layer, the first and second adjacent layers, and the first, second, and third adjacent layers as the range. These two types of diffusion circles are used as both a microenvironment consideration and important local information for calculating the interaction index of the muscle fiber.

[0145] S3: Connected domain microenvironment extraction: Based on the adjacency matrix of muscle fibers, the connected domains formed by the adjacency relationship of the same type of muscle fibers are identified, which serves as both a microenvironment consideration and an important local information for calculating the interaction index of the muscle fiber.

[0146] S4: Muscle fiber interaction feature measurement: measures the degree of contact between a muscle fiber and its various types of muscle fibers, measures its multi-layered spatial interaction information in the diffusion circle, and measures its centrality in the connected domain as important spatial interaction features.

[0147] S5: Microenvironment characteristic indicators: Measures the size and shape of connected domains, the composition ratio and distribution pattern of various types of muscle fibers in the diffusion circle, and statistical measures of muscle fiber properties in these two microenvironments.

[0148] S6: Microenvironment classification: Similarity clustering is performed on diffusion circles and connected regions based on rotational overlap, weighted Euclidean distance of microenvironment feature vectors, and weighted k-means and PCA k-means clustering methods.

[0149] S7: Measurement of the changing trend of diffusion circle characteristics with increasing radius: For each muscle fiber, a diffusion circle is formed at 1.5, 2.5, and 3.5 times the average diameter of all muscle fibers, and the first adjacent layer, the first and second adjacent layers, and the first, second, and third adjacent layers are used as the range. The proportional vector index in the three adjacent diffusion circles is generated and the first and second differences are calculated.

[0150] The steps for constructing the distance matrix in S1 include the following operations:

[0151] Calculate the center distance and shortest distance between any two muscle fibers in the image, where the center distance between the i-th muscle fiber and the j-th muscle fiber is denoted by cd. ij express:

[0152]

[0153] Shortest distance is MD ij express:

[0154]

[0155] Where Ω i and Ω j Let represent the sets of contour coordinates of the i-th and j-th muscle fibers, respectively. The distance matrix CD = cd ij Distance matrix MD = {md ij}

[0156] Extraction of the diffusion circle microenvironment in S2 includes the following operations:

[0157] S2.1: Construct an adjacency matrix for all muscle fibers. Initialize the adjacency matrix ADJ to the distance matrix MD. For elements whose values ​​are less than 1 / 2 of the minimum length of the minor axis of all muscle fibers (denoted as MSD), set their values ​​to 1, thus completing the construction of the adjacency matrix ADJ.

[0158] The adjacency matrix ADJ has elements that satisfy:

[0159]

[0160] Where md ij The shortest distance matrix for muscle fibers is MD = {md} ij The element in} represents the shortest distance between the i-th and j-th muscle fiber contours, and MSD is the minimum minor axis of all muscle fibers in the image.

[0161] S2.2: The first type of diffusion circle is uniquely represented (determined) by the center point of the microenvironment and the absolute radius of the microenvironment. The center point of the microenvironment is the geometric center of a certain muscle fiber, which is matched with the muscle fiber ID. At the same time, circles whose contours exceed the image boundary are excluded. For each valid diffusion circle with three radius lengths, the muscle fiber ID where the center is located, the corresponding coordinates, and the radius of the diffusion circle are extracted as the identifier of this microenvironment.

[0162] S2.3: The second type of diffusion circle is uniquely represented (determined) by the center point of the microenvironment and the adjacency radius of the microenvironment. The adjacency radius is defined by the distance between muscle fibers in the adjacency matrix. If muscle fiber i needs to move k steps to reach muscle fiber j through adjacency relationships, then their distance is k. For each valid diffusion circle (excluding circles whose contours exceed the image boundary) with three radius lengths of adjacency distances of 1, 2, and 3, the muscle fiber ID at the center, the coordinates of the center, and the adjacency radius of the diffusion circle are obtained as the identifier of this microenvironment. The muscle fiber IDs contained within it are calculated. Taking the diffusion circle of the i-th muscle fiber with an adjacency radius of k as an example, the calculation method for the muscle fiber IDs contained within it is as follows:

[0163] find(ADJ(i,:)+ADJ 2 (i,:)+…+ADJ k (i,:))

[0164] ADJ k (i,:) represents the i-th row of the adjacency matrix raised to the power of k, and the find function returns the position of the non-zero element in it.

[0165] The extraction of connected component microenvironments in S3 includes the following operations:

[0166] S3.1: Construct the adjacency matrix ADJ for all muscle fibers of class k. (k) ADJ (k) =ADJ(Nk N k ), where N k Let ID (set of numbers) be the ID of all k-th type muscle fibers in the image.

[0167] S3.2: Connected Component Identification: For each type of muscle fiber (taking the k-th type as an example), the breadth-first search algorithm is used to identify connected components.

[0168] Initialization: Create a visited array to record visited nodes, and create a connected_components list;

[0169] Traversal: Traverse the adjacency matrix ADJ (k) For each row i, for each unvisited i, perform a BFS (Breadth-First Search) traversal and create a new connected set;

[0170] BFS Traversal: Define an empty queue, add the starting node `i` to the queue, and add the current node `i` to the head of the current connected vector (storing the indices of nodes in the connected set). When the queue is not empty, remove the first node `current_node` from the queue, mark it as visited, and add it to the current connected set. Then traverse all adjacent nodes of `current_node`. If a neighbor has not been visited and is adjacent to the current node, add the neighbor node to the queue for later access, mark it as visited, and add it to the current connected vector. The loop continues until the queue is empty. Whenever a connected set is found, add it to the `connected_components` list. After traversal, sort each set in the list by the number of elements in descending order and return the list. The `i`th element in the list represents the set of myofibrous IDs in the `i`th connected component. The specific process is as follows:

[0171]

[0172] visited[i] = True / / Mark node i as visited

[0173] While queue is not empty: / / Perform BFS traversal while queue is not empty.

[0174] current_node = Dequeue(queue) / / Retrieve the first node from the queue

[0175] For j from 0 to size-1: / / Iterate through all neighboring nodes of current_node (j is the neighbor node ID).

[0176] If visited[j] is False and adj_matrix[current_node][j] is 1: / / Neighbors have not been visited and are connected by edges.

[0177] Enqueue(queue, j) / / Adds neighboring nodes to the queue for later access.

[0178] visited[j] = True / / Mark neighbor nodes as visited

[0179] current_connected_component.add(j) / / Add the neighbor node to the current connected vector.

[0180] connected_components.add(current_connected_component) / / Add the connected components to the result list after the connected components are generated.

[0181] Sort connected_components in descending order of the number of elements in the set / / Sort by the size of the connected set from highest to lowest

[0182] Return connected_components / / Returns a sorted list of connected components

[0183] End Function

[0184] The measurement of myofibril interaction characteristics in S4 includes the following operations:

[0185] S4.1: Calculate the ratio of the contact circumference of each muscle fiber to that of each type of muscle fiber, including the following steps:

[0186] S4.1.1: Filtering adjacent muscle fibers (taking the i-th fiber as an example): Extract the i-th row of ADJ (excluding the current muscle fiber i itself), and extract the set of k-th type muscle fiber IDs in the i-th row of ADJ that are adjacent to the target muscle fiber i according to the muscle fiber type.

[0187] S4.1.2: Set a threshold for the distance between adjacent pixels: for example, one-third of the minimum shortest diameter of the global muscle fiber;

[0188] S4.1.3: For each category (taking the k-th category as an example), traverse adjacent muscle fibers to calculate contact pixels: For each muscle fiber j adjacent to the i-th muscle fiber, first count the number of pixels whose shortest distance from all pixels on the i-th contour to pixels on the j-th contour is less than the above threshold, denoted as n. ij Count the number of pixels on contour j whose shortest distance to a pixel on contour i is less than the above threshold, and denote this as n. jiThen, n ij Divide this by the total number of pixels in contour i to obtain the "adjacent pixel ratio" of i to j; then divide n... ji Divide by the total number of pixels in the j-contact profile to obtain the "adjacent pixel ratio" of j to i;

[0189] S4.1.4: Summing to obtain the final ratio: Add the "outgoing ratio of adjacent pixels" of all adjacent muscle fibers j to get the "outgoing ratio" of muscle fiber i; add the "ingoing ratio of adjacent pixels" of all adjacent muscle fibers j to get the "ingoing ratio" of muscle fiber i. The "outgoing ratio" and "ingoing ratio" describe the channel width of the interaction between muscle fibers and various types of muscle fibers;

[0190] The calculation formula is as follows:

[0191] The ratio of out-of-pixel distances between muscle fiber i and its neighboring k-th type muscle fiber, and the ratio of in-of-pixel distances between the k-th type muscle fiber and its neighboring i, are respectively:

[0192]

[0193] Where N(i) (k) It is the set of the kth type of muscle fibers adjacent to muscle fiber i, |c i | and |c j | represents the total number of pixels in the contour of muscle fibers i and j (reflecting the perimeter scale).

[0194] S4.2: Calculate the weighted sum of the contact ratios between each muscle fiber and surrounding muscle fibers, selecting specific weights to reflect the interaction strength at a specific level, through:

[0195]

[0196] Where N(i) is the set of muscle fibers adjacent to muscle fiber i, W j It is a specified attribute value of the j-th muscle fiber, such as area, expression level of a gene, or activity of a pathway.

[0197] S4.3: Calculate the information within the circle radiating outward from the muscle fiber (divided into absolute diffusion circles and adjacent diffusion circles; the calculation methods for the following features are basically the same for both types of diffusion circles). Taking the k-th muscle fiber as an example, the information of the diffusion circle centered on it is condensed into a matrix consisting of distance, included angle, category (membership degree), and occupied area, calculated as follows:

[0198] S4.3.1: For an absolute diffusion circle, if the geometric center of a muscle fiber is inside the circle, it is included in the statistics; for adjacent diffusion circles, the muscle fibers contained therein are taken.

[0199] S4.3.2: Extract the center point of these muscle fibers (taking the k-th muscle fiber as an example, using a diffusion circle with radius (number of adjacent layers) as r). Category (Membership Degree) i represents the i-th muscle fiber within the diffusion circle;

[0200] S4.3.3: Calculate the area of ​​intersection between these muscle fibers and the diffusion circle:

[0201]

[0202] in It represents the area of ​​the diffusion circle, and r represents the radius of the diffusion circle or the number of adjacent concentric rings. This represents the area of ​​the i-th muscle fiber;

[0203] S4.3.4: Calculate the distance and angle from the center point of these muscle fibers to the center of the circle (excluding the muscle fibers corresponding to the center of the circle);

[0204] distance:

[0205] Included angle: First calculate the radians:

[0206] Among them O (k,r) =(x o ,y o () represents the coordinates of the center of the diffusion circle. Let sign(y) represent the geometric center coordinates of the i-th muscle fiber within the diffusion circle. i -y o ) is a symbolic function:

[0207]

[0208] Then calculate the included angle:

[0209] S4.3.5: Concatenate columns to output the adjacency information matrix.

[0210] The adjacency information matrix describes the adjacency of the muscle fiber, while the adjacency with gradually increasing radius describes the change in the distribution of surrounding muscle fibers as the distance increases.

[0211] S4.4: Measure the interaction strength between muscle fibers and different types of muscle fibers in a specified local area (within a specified radius or adjacent concentric layers):

[0212] Taking the i-th muscle fiber as an example, the interaction strength between it and any j-th muscle fiber within a diffusion circle of a specified radius or adjacent layer r is measured using a gravitational model as follows:

[0213]

[0214] Where A i Let CD be the area of ​​the i-th muscle fiber (or it can be replaced by the class membership [-1, 1], where -1 represents type II muscle fiber and 1 represents type I muscle fiber, as well as the activity / expression level of certain genes or pathways), CD = cd ij Let be the center distance between the i-th and j-th muscle fibers. The assumption behind this setting is that the larger the area, the richer the information it contains, and the stronger its ability to receive or transmit information. The farther the distance, the more difficult it is for the muscle fibers to communicate signals and molecules with each other. One reason for using the center distance instead of the shortest distance is that the shortest distance cannot avoid the influence of shape, such as for long and thin muscle fibers, while the center distance can better improve this.

[0215] The sum of the interaction strengths between the i-th muscle fiber and any k-th type of muscle fiber within the diffusion circle is:

[0216]

[0217] in Represents the set of IDs of the k-th type of muscle fiber within a diffusion circle with radius r or adjacent concentric circles of the i-th muscle fiber;

[0218] S4.5: Measures the centrality of a muscle fiber's position within a connected region formed by similar muscle fibers.

[0219] For each type of muscle fiber (taking the k-th type of muscle fiber as an example), extract the adjacency matrix of the connected components from the total adjacency matrix Adj. (Taking the i-th connected region of the k-th type of muscle fiber as an example), then output the proportion of adjacent pixels between muscle fibers within the connected region. Opposite edge (N) i (k) (j)→N i (k) (j) is weighted according to the area of ​​the muscle fiber. For node N i (k) (j) Perform weighting to form a weighted adjacency matrix.

[0220]

[0221] Where N i (k) This represents the set (vector) of muscle fiber IDs in the i-th connected region of the k-th type of muscle fiber; the other symbols are defined above.

[0222] Using a weighted adjacency matrix Calculate the weighted closeness centrality, weighted betweenness centrality, weighted degree centrality, and weighted eigenvector centrality for each muscle fiber in the connected component. For the v-th muscle fiber in the connected component, its four centralities are calculated as follows:

[0223] The weighted closeness centrality is calculated as follows:

[0224]

[0225] Where d(u,v) is the shortest weighted distance between muscle fibers u and v in the weighted adjacency matrix, V represents the set of all muscle fibers in the connected region, and n is the number of muscle fibers in the connected region. Weighted closeness centrality measures how close a given muscle fiber in a connected region is to all other muscle fibers in the region. A muscle fiber with high weighted closeness centrality means it can reach other muscle fibers in the connected region via a shorter path, potentially offering advantages in signal, molecular propagation, and reception.

[0226] Weighted betweeness centrality:

[0227]

[0228] Where σ st σ is the total number of shortest paths between muscle fibers s and t in the weighted adjacency matrix. st (v) represents the number of shortest weighted paths between muscle fibers s and t that pass through muscle fiber v. Betweenness centrality measures the frequency with which a given muscle fiber appears on the shortest paths between other pairs of muscle fibers within a connected region. Muscle fibers with high weighted betweenness centrality may play a crucial bridging role in information transmission, resource transfer, and other processes.

[0229] Weighted in-degree centrality:

[0230]

[0231] Where k v is the degree of muscle fiber v (i.e., the number of muscle fibers directly connected to muscle fiber v in the adjacency matrix), and n is the total number of muscle fibers in the connected domain.

[0232] Weighted in-degree centrality is a simple metric for measuring the importance of muscle fibers in a connected domain. It refers to the sum of the edge weights pointing to that muscle fiber. A high in-degree centrality indicates that the muscle fiber has a greater degree of direct spatial contact with other muscle fibers of the same type.

[0233] Weighted eigenvector centrality:

[0234]

[0235] Here, x is the eigenvector, and the L1 norm (sum of absolute values) of x is constrained to 1, i.e., |x|1 = 1. This condition is used to eliminate scale ambiguity in the eigenvectors and ensure the comparability of centrality values. λ is the largest eigenvalue. The element values ​​of the eigenvector x corresponding to the largest eigenvalue are the eigenvector centralities of each node. The weighted eigenvector centrality considers the importance of the muscle fibers' neighboring muscle fibers within the connected domain. The importance of a muscle fiber depends not only on the number of muscle fibers it is directly connected to but also on the importance of these neighboring muscle fibers. Its value reflects the frequency with which information passes through that muscle fiber (considering only spatial factors).

[0236] Note that this solution can also provide an unweighted version of the above-mentioned muscle fiber centrality metric, that is, without weighting the edges and nodes of the adjacency matrix, with the same calculation method.

[0237] The S5 microenvironment characteristic index measurement includes the following operations:

[0238] S5.1: Quantification of diffusion circle microenvironment characteristics, including:

[0239] S5.1.1: The proportion of the total area of ​​each type of muscle fiber within the diffusion circle (proportion vector). If the obtained value is a category membership degree, it is classified according to a threshold of 0.5 (greater than or equal to 0.5 represents type I, less than or equal to 0.5 represents type II). Taking the i-th muscle fiber and the diffusion circle with radius or adjacent concentric circles of r as an example, it is denoted as ARatio. i r .

[0240] S5.1.2: The proportion of arc occupied by each type of muscle fiber within the diffusion circle (proportional vector), taking the k-th type as an example (the muscle fiber located at the center of the diffusion circle should be excluded when calculating the arc), includes the following steps;

[0241] S5.1.2.1: First, define the center O of the current diffusion circle. (k) =(x o ,y o The set of contour points of the k-th type of muscle fiber: n is the number of contour points, n≥0.

[0242] S5.1.2.2: Polar Angle Calculation and Standardization: For each contour point Calculate relative to the center O (k) polar angle And standardize to the interval [0, 2π):

[0243]

[0244] S5.1.2.3: Calculation of maximum sector radian: If n≤1 (no valid sector or single point), the radian is 0; otherwise, calculate the maximum value of the absolute values ​​of the differences between all pairs of polar angles:

[0245]

[0246] S5.1.2.4: Calculation of Radius Ratio: The radius ratio of the k-th type of muscle fiber within the diffusion circle is...

[0247] This radian scale vector This reflects the degree of dispersion of different types of muscle fibers within the diffusion circle (taking the i-th muscle fiber as an example, with a radius or adjacent concentric circles of r). A smaller value indicates a greater concentration within a small area of ​​the myocardial fibers, while a larger value indicates a wider distribution. For example, Figure 1 (A) In the diffusion circle marked in the sub-icon, type I muscle fibers occupy a large proportion of the arc, exceeding 90%. Figure 1 (B) The proportion of type I muscle fiber curvature in the diffusion circle marked by the sub-icon is relatively small, less than 70%;

[0248] S5.1.3: The adjacency information matrix of the i-th muscle fiber radius or adjacency layer r constructed above.

[0249] S5.1.4: Statistical values ​​of interaction attributes (mean, variance, median) for different types of muscle fibers;

[0250] S5.2: Calculate the autocorrelation Moran coefficient of the muscle fiber index values ​​(using area index as an example) within the diffusion circle. This step calculates the values ​​for all muscle fibers within the diffusion circle and for only a specific type of muscle fiber (e.g., type I muscle fibers). Specifically, it includes the following steps:

[0251] S5.2.1: Based on the current muscle fiber ID, obtain its geometric center coordinates. Using this geometric center as the center, and combining it with a given absolute radius or adjacent radius, determine all muscle fiber IDs and their corresponding contour coordinates within the diffusion circle. To focus on a specific type of muscle fiber, further filter out the muscle fiber IDs and their contour coordinates belonging to that type from the muscle fibers within the aforementioned diffusion circle.

[0252] S5.2.2: Extract the center distance matrix CD{cd based on these muscle fiber IDs. ij Information;

[0253] S5.2.3: The inverse distance weighting method is adopted. When the distance is 0, the weight is set to 0 (excluding itself). Otherwise, the weight is the reciprocal of the distance (the closer the distance, the higher the weight).

[0254] Mathematical expression:

[0255] cd ij w is the center distance between muscle fibers i and j. ij As weight;

[0256] S5.2.4: Center the target feature values ​​(such as area, class membership) of muscle fibers within the circular region:

[0257]

[0258] Where n is the total number of muscle fibers in the world. Let x be the feature value of the i-th muscle fiber in the global range. i Let z be the characteristic value of the i-th muscle fiber within the diffusion circle. i These are the local eigenvalues ​​after centering;

[0259] S5.2.5: Calculation of local Moran coefficients:

[0260]

[0261] Where m is the number of muscle fibers within the diffusion circle;

[0262] S5.2.6: Calculate the P-value for the substitution test:

[0263]

[0264] By randomly permuting the eigenvalues ​​N (set to 999) times, the distribution of the Moran coefficient under the null hypothesis was simulated, and the significance of the observed values ​​was calculated. The Moran index characterizes the spatial autocorrelation of myofiber indicators within the microenvironment, and more phenotypic indicators can be added later to discover the promoting and inhibiting effects between myofibers;

[0265] S5.2.7: Local Moran Coefficient and Global Moran Coefficient (I global (Compare)

[0266] To determine whether the local Moran coefficient exhibits clustering (positive correlation) or exclusion (negative correlation) relative to the global Moran coefficient, it is necessary to combine statistical significance with the actual value of the local Moran coefficient.

[0267]

[0268]

[0269] S5.3: Calculate the local distance matrix entropy and angle matrix entropy values ​​for each type of muscle fiber (taking type I as an example) within the diffusion circle, including:

[0270] S5.3.1: First, extract the ID, geometric center coordinates, and contour coordinates of the type I muscle fibers inside the effective diffusion circle.

[0271] S5.3.2: Extract the local center distance matrix CD based on muscle fiber ID. (1) Calculate the included angle matrix A (1) The included angle matrix is ​​calculated as follows:

[0272] S5.3.2.1: The input geometric center coordinate matrix is ​​C∈R n×d Where n is the number of geometric centers, and d is the coordinate dimension of each geometric center. The i-th row of C... i The coordinate vector representing the i-th geometric center, i.e.

[0273] S5.3.2.2: Calculate the dot product matrix: Dot product matrix D∈R n×n Element D in ij c represents the vector of the i-th geometric center. i and the j-th geometric center vector c j The dot product is calculated using the following formula:

[0274] S5.3.2.3: Calculate the L2 norm of each geometric center: for each geometric center vector c i Its L2 norm (also known as the Euclidean norm) is defined as:

[0275]

[0276] S5.3.2.4: Constructing the included angle matrix: Included angle matrix A (1) ∈R n×n The element in c represents the vector of the i-th geometric center. i and the j-th geometric center vector c j The cosine of the included angle, according to the definition of the vector dot product:

[0277]

[0278] S5.3.2.5: Since the angle between each vector and itself is 0 degrees, and the cosine value is 1, it is necessary to adjust the angle matrix A. (1) The diagonal element is set to 1, i.e., A ii =1, i = 1, 2, ..., n;

[0279] S5.3.3: Convert the distance matrix and the included angle matrix into vectors and remove the zero values ​​to obtain the distance matrix vector d and the included angle matrix vector a.

[0280] S5.3.4: Data Discretization: Divide vectors d and a into k intervals, with k initially set to 20. Let the intervals be B1, B2, B3, ..., B kThese intervals are generated using the `np.linspace` function. `bins = np.linspace(min(v), max(v), k+1)`, where `v` can be either `d` or `a`.

[0281] S5.3.5: Calculate the frequency of intervals: For each vector, count the number of elements within each interval. Let f i For the i-th interval B i The frequency of the inner element, i = 1, 2, ..., k.

[0282] S5.3.6: Calculate the probability distribution: Remove intervals with a frequency of 0, then calculate the probability p for each non-empty interval. i :

[0283]

[0284] S5.3.7: Calculate the entropy value: This allows us to obtain the entropy value of the distance matrix between type I muscle fibers within the diffusion circle. Entropy of the included angle matrix

[0285] The meaning of entropy: The larger the value, the more scattered the distance distribution between type I muscle fibers within the diffusion circle, with different distances appearing at similar frequencies. The smaller the value, the more concentrated the distance between pairs of type I muscle fibers within the diffusion circle is in a few intervals, indicating a lower degree of disorder. The larger the value, the more uniform the directional angle of the geometric center vector of the type I muscle fibers within the diffusion circle is, indicating that the orientation of each muscle fiber has no obvious dominant direction and the arrangement is disordered. The smaller the value, the more similar the geometric center vectors of most muscle fibers in type I within the diffusion circle are, meaning that the muscle fibers tend to be arranged in the same direction and have a strong orientation consistency.

[0286] S5.4: Measure the degree of aggregation of different types of muscle fibers within the diffusion circle (taking the k-th type of muscle fiber in the diffusion circle with radius r of the i-th muscle fiber or adjacent concentric circles as an example).

[0287] S5.4.1: Extracting features of all types of muscle fibers within the diffusion circle: Identify diffusion circle C i All muscle fibers are included, and the geometric center coordinates p of each muscle fiber are extracted. j Calculate the area A of each muscle fiber. j Extract the geometric center coordinates p of each muscle fiber of class k from it. j And the area A of each muscle fiber j .

[0288] S5.4.2: Calculate the area-normalized K-type distance ratio (AKDR):

[0289] First: Calculate the normalized average distance of the area within the k-th type of muscle fiber:

[0290] Where d(p) j ,p l N represents the Euclidean distance between two geometric centers. k This represents the number of the k-th type of muscle fibers within the diffusion circle. The reason for using area as the denominator for normalization is to eliminate the distance effect caused by different areas;

[0291] Then, calculate the area-normalized mean distance of all muscle fibers:

[0292] Where N is the number of all muscle fibers within the diffusion circle;

[0293] S5.4.3: Calculate the final degree of aggregation: when When, it indicates that the aggregation degree of the k-th type of muscle fiber within the diffusion circle is relatively higher; when When the value is 0, it indicates that the aggregation degree of the k-th type of muscle fiber within the diffusion circle is relatively lower, and it tends to be distributed in a dispersed manner.

[0294] S5.5: Quantification of microenvironment features of connected domains, including the following steps:

[0295] S5.5.1: Quantitative Size Indicators: Area and Number of Muscle Fibers

[0296] S5.5.1.1: Extract the IDs of all muscle fibers in the connected component, and then calculate the area of ​​each muscle fiber. For the area a of the v-th muscle fiber in the connected component... v :

[0297]

[0298] Where n represents the number of coordinate points of the muscle fiber contour.

[0299] S5.5.1.2: Summing the areas of the muscle fibers within each connected region yields the area A of the connected region. total : Where m represents the number of muscle fibers within the connected domain.

[0300] S5.5.1.3: Finally, calculating the length of the list of muscle fibers contained in the connected region yields the number of muscle fibers. The area of ​​the connected region and the number of muscle fibers reflect the size of the region it occupies.

[0301] S5.5.2: Quantitative Shape Indicators: Aspect Ratio, Roundness, Fillness (Firmness), Rectangularity

[0302] Before calculating the shape index of the connected components, the contours of the muscle fibers within the connected components are first merged, utilizing the uniqueness of the convex hull: for a given set of points, its convex hull is a unique minimal convex polygon. The process is as follows:

[0303] S5.5.2.1: Extract the set of contour coordinates {Ω1,Ω2,…,Ω} for each muscle fiber within the connected region. m},

[0304] S5.5.2.2: Merge contour coordinate point sets

[0305] S5.5.2.3: Calculate the vertices of the convex hull: Vertices = ConvHull(Ω) merged ),

[0306] S5.5.2.4: Output: Merged contours (convex hull vertex sequence) Vertices.

[0307] This process ensures that the merged contour is the smallest convex polygon containing all the original contour points, providing a unified outer boundary representation for subsequent calculations of connected domain geometric features.

[0308] Aspect ratio calculation: Obtain the width w and height h using cv2.minAreaRect(Vertices). A larger value indicates that the connected region has a more elongated shape.

[0309] Circularity: Where L is the perimeter of the connected region, obtained using cv2.arcLength(Vertices, True). The closer the value is to 1, the closer the shape of the connected region is to a circle.

[0310] Fill factor: is the sum of the areas of muscle fibers within a connected region divided by the area of ​​the convex hull of the merged contour of the muscle fibers within the connected region. Where A Vertices It is obtained from cv2.contourArea(contour). The closer the fill value is to 1, the more compact the distribution of muscle fibers within the connected region.

[0311] Rectangularity: The closer the value is to 1, the closer the shape of the connected component is to a rectangle.

[0312] S5.5.3: Cluster Coefficient

[0313] First, extract or calculate the adjacency matrix ADJ of the connected components. c Through ADJc Construct an undirected graph G = (V, E), where the vertex set V = C = {v1, v2, ... v} n}, where n is the total number of muscle fibers in the connected region, and the edge set E is determined by ADJ. c The zero element in the middle and non-zero elements has been determined.

[0314] (1) For any node v in graph G, its clustering coefficient is:

[0315] Where, k v E is the degree of node v. v : The actual number of edges between all adjacent nodes of node v (i.e., the number of edges between adjacent node pairs in the graph).

[0316] (2) The clustering coefficient of the connected components is the average of the clustering coefficients of all nodes in the graph:

[0317]

[0318] The clustering coefficient measures the distribution of muscle fibers within a connected region; a larger value indicates a more clustered distribution. A schematic diagram illustrating the diffusion circle, connected regions of similar muscle fibers, and some basic features is shown below. Figure 2 .

[0319] S6 microenvironment classification includes the following operations:

[0320] S6.1: For absolute diffusion circles, a sparse matrix clustering method is used to reduce storage and time consumption. Sparse matrix clustering includes the following steps:

[0321] S6.1.1: Constructing the diffusion circle category value matrix: For all effective circular regions generated by muscle fibers in each image, the muscle fibers inside are truncated based on their membership degree [0,1]. For example, membership degrees greater than or equal to 0.5 are type I, and those less than 0.5 are type II (this can be changed according to the actual situation). The diffusion circle... Pixels within the contour range of type I muscle fibers are set to 1, pixels within the contour range of type II muscle fibers are set to 2, and all other pixels are set to 0.

[0322]

[0323] Where i represents the i-th muscle fiber (i = 1, 2, ..., n is the total number of effective muscle fibers in the image), r represents the radius of the diffusion circle, and (x, y) are the coordinates of the pixels within the rectangular region corresponding to the circular domain. The circumscribed rectangle of this diffusion circle is then used to generate a matrix M of the same size. i , and initialize it to 0, then let

[0324] S6.1.2: Initialize the similarity matrix: Create a circular domain similarity matrix and initialize all elements in the matrix to 0. For an n×n similarity matrix S (where n is the number of effective diffusion circles), we have: S ij =0;

[0325] S6.1.3: Calculate the circular region similarity matrix: Compare the pixel values ​​of the circular regions pairwise. To find the best match, one of the diffusion circles is rotated during the comparison process. Specifically, the first diffusion circle is fixed, and the second diffusion circle is rotated by 0 degrees, 90 degrees, 180 degrees, and 270 degrees respectively. Each rotation is then compared with the first diffusion circle, and the number of identical pixel values ​​is counted. The maximum value among these comparison results is taken as the number of identical pixels between the two circular regions, and this value is written into the corresponding position in the similarity matrix.

[0326] First, define the function g(M) a M b This is used to calculate the number of corresponding pixel values ​​that are equal in two matrices.

[0327] g(M a M b )=Σ x Σ y δ(M a (x,y),M b (x,y))-N O ,

[0328] Where N represents O The number of pixels located outside the diffusion circle in this matrix, δ(a, b) is the Kronecker function:

[0329]

[0330] M j Rotating by 0°, 90°, 180°, and 270° respectively yields... and Then M i and M j Similarity S ij for:

[0331]

[0332] Where, N circle This represents the number of pixels within the diffusion circle.

[0333] S6.1.4: K-means is based on the logic that "the smaller the distance, the more similar the samples". After calculating the pairwise similarity matrix, it is converted into a distance matrix. Let the distance matrix be D∈R. n×n The element is defined as:

[0334] D[i,j]=1-S ij ,

[0335] Among them, S ij These are the elements at the corresponding positions in the original similarity matrix.

[0336] S6.1.5: Multidimensional Scaling (MDS) Transformation:

[0337] The distance matrix D is mapped to a low-dimensional feature matrix X∈R using the MDS algorithm. n×d (where d is a low-dimensional dimension), satisfying the condition that for any i, j, X, the Euclidean distance between samples i and j is approximately equal to the original distance:

[0338] ||X[i,:]-X[j,:]||2≈D[i,j]

[0339] In the formula, X[i,:] represents the i-th row of the feature matrix (the low-dimensional eigenvector of the i-th circle);

[0340] S6.1.6 The elbow rule determines the number of clusters, optimal_clusters, and the clustering results. The maximum number of clusters, num_clusters, is set to 20.

[0341] For each cluster size k cluster Solve for the objective function:

[0342]

[0343] C i It is the i-th cluster, u i It is cluster C i The geometric center (obtained by calculating the X matrix).

[0344] Calculate the second derivative of SSE:

[0345] Where C i For the i-th cluster, u i It serves as the cluster center.

[0346] SSE″(k)≈SSE(k+1)-2SSE(k)+SSE(k-1),

[0347] Where k∈[2,num_clusters-1], k-1=0 cannot be calculated when k=1; k+1=num_clusters cannot be calculated when k=num_clusters (it exceeds the initially set range).

[0348] Optimal number of clusters:

[0349]

[0350] S6.1.7 Perform K-means clustering for a given number of clusters.

[0351] Assuming the number of clusters is k, weighted K-means clustering is performed on the low-dimensional feature matrix X, as follows:

[0352] 1. Initialization: Randomly select k samples as initial cluster centers. (Each center is a d-dimensional feature vector);

[0353] 2. Iterative optimization:

[0354] Sample allocation: For each row X[i,:] of X, calculate its weighted distance to the cluster centers and allocate it to the cluster with the smallest distance: Where (t is the iteration number, Let i be the sample set of the i-th cluster at the t-th iteration. (It is the coordinate of its cluster center).

[0355] Update cluster centers: Recalculate the cluster centers (weights affect mean calculation): ( (where is the number of samples in the i-th cluster at the t-th iteration).

[0356] Convergence condition: Iteration stops when the change in cluster centers is less than a given threshold ε.

[0357]

[0358] Output: Cluster label l for each circle i ∈{1,2,…,k}.

[0359] S6.1.8 Calculation of average similarity within clusters (quantifying cluster quality)

[0360] Based on the original normalized similarity matrix, the average similarity within each cluster is calculated to verify that "samples of the same type have high similarity":

[0361] 1. Extract the local similarity matrix: For the C-th cluster (C = 1, 2, ..., ooptimal_clusters), let its contained circle index be... (n C Extract a submatrix from the similarity matrix (where the number of samples within the cluster is 1).

[0362] 2. Calculate the upper triangular average: Since local_S_c is symmetric and the main diagonal (i=j) is meaningless, the upper triangular portion is used to calculate the average: (mean_sim_c∈[0,1], the larger the value, the higher the similarity of the circles within the cluster, and the better the clustering effect).

[0363] S6.2: Use weighted k-means and PCA-Kmeans clustering for connected components.

[0364] Weighted k-means:

[0365] (1) Data Reading and Feature Normalization: The pandas library is used to read the data from the document storing connected component information. For the features that need to be clustered, features = ['circularity', 'aspect ratio', 'fill degree', 'rectangularity', 'total area', 'cluster coefficient'], MinMaxScaler is used for normalization. The normalization formula is:

[0366]

[0367] Here, x is the original feature value, that is, the specific numerical value of a sample in the dataset for a certain feature. For example, for the feature "circularity", the circularity value of a connected component sample is x. min and x max These are the minimum and maximum values ​​of the feature, x. norm These are the normalized eigenvalues.

[0368] (2) Weight settings: Custom weight list w * = [0.1, 0.1, 0.1, 0.1, 0.5, 0.1]. Since the total area is more intuitive, the weight is set to the maximum (this can be changed according to specific needs).

[0369] (3) The elbow rule determines the number of clusters optimal_clusters and the clustering results, and the maximum number of clusters num_clusters is set to 20.

[0370] For each k cluster Solve for the objective function:

[0371]

[0372] Where, c1,…, : Represents different clustering partitioning schemes, with the goal of finding the optimal partitioning scheme that minimizes the objective function value. Double summation: outermost layer Summing over all clusters, intermediate layer Sum the squared distances between all samples within each cluster and the cluster center. Weighted Euclidean distance:

[0373]

[0374] in, Let u be the normalized eigenvector of the j-th connected component. iLet be the center (mean) of the eigenvectors of the connected components in the i-th cluster. The weight of the k-th feature.

[0375] The objective function represents the minimum sum of the weighted squared distances from all samples to the center of their respective clusters under different clustering schemes.

[0376] Calculate the second derivative:

[0377]

[0378] SSE″(k)≈SSE(k+1)-2SSE(k)+SSE(k-1),

[0379] Where k∈[2,num_clusters-1], k-1=0 cannot be calculated when k=1; k+1=num_clusters cannot be calculated when k=num_clusters (it exceeds the initially set range).

[0380] Optimal number of clusters:

[0381]

[0382] PCA-Kmeans clustering:

[0383] (1) Data Reading and Feature Selection: Taking connected components as an example, read the connected component data and select 6 original features for analysis: features = ['circularity', 'aspect ratio', 'fill degree', 'rectangularity', 'total area', 'clustering coefficient'], and construct the original feature matrix X. Let the number of connected components of an image be n, the number of features be p = 6, and the original feature matrix X = [x ij ]∈R n×p , where x ij It represents the j-th eigenvalue of the i-th connected component.

[0384] (2) Z-score standardization (data preprocessing): Standardize the original features to eliminate the influence of dimensions and make each feature follow a standard normal distribution with a mean of 0 and a standard deviation of 1.

[0385] in

[0386] (3) Principal Component Analysis (PCA dimensionality reduction):

[0387] First, calculate the eigenvalues ​​and eigenvectors of the covariance matrix of the standardized data.

[0388] Covariance matrix:

[0389] Calculation of eigenvalues and eigenvectors: Solve ∑v k = λ k v k , where λ k is the eigenvalue (in descending order), and v k is the corresponding unit eigenvector (2-norm is 1).

[0390] Then, calculate the ratio of explained variance. Proportion of explained variance by a single principal component:

[0391]

[0392] where λ k is the eigenvalue corresponding to the k-th principal component (the k-th eigenvalue), p is the total number of original features, sum of all principal component eigenvalues.

[0393] Cumulative proportion of explained variance:

[0394] where k represents the first k principal components (k ≤ p);

[0395] The cumulative proportion of explained variance represents the proportion of the total variance of the original data jointly explained by the first k principal components;

[0396] Next, determine the number of principal components. When n ≥ p, select k such that the cumulative variance ratio ≥ 85%; when n < p, select k = min(n - 1, p) and the cumulative variance ratio ≥ 80% (can be adjusted flexibly. If the cumulative variance ratio of the first k principal components in the actual data is close to the threshold, but the improvement is small after adding 1 principal component).

[0397] Calculate the principal component scores. Project the standardized data onto the principal component space to obtain the principal component score matrix Z after dimensionality reduction: Z = X<00​​​​​​​​​​​​​​​​​​​​​​​​​​​To determine the clustering scheme, where C i Let u be the set of indices contained in the i-th cluster. i Let Z be the center vector of Z in the Z matrix, which is the mean of the corresponding principal component score vector.

[0403] Calculate the second derivative:

[0404]

[0405] SSE″(k)≈SSE(k+1)-2SSE(k)+SSE(k-1),

[0406] Where k∈[2,num_clusters-1], k-1=0 cannot be calculated when k=1; k+1=num_clusters cannot be calculated when k=num_clusters (it exceeds the initially set range).

[0407] Optimal number of clusters:

[0408]

[0409] S6.3: For valid adjacent diffusion circles, the above PCA-Kmeans clustering is used. The selected features are [the proportion of the total area of ​​each type of muscle fiber (proportion vector), the proportion of the arc occupied by each type of cell (proportion vector), the entropy values ​​of the local distance matrix and the included angle matrix of each type of muscle fiber, and the degree of clustering of each type of muscle fiber], which are then stretched into a feature vector instead of a matrix. If the adjacency radius exceeds 1, for example, k, then the features of the adjacent diffusion circle are calculated for each adjacency radius i = 1, 2, ..., k, and these features are then concatenated into a feature vector. Next, based on these feature vectors, the above PCA-Kmeans clustering is performed.

[0410] The measurement of the change trend of the diffusion circle feature with increasing radius or adjacent concentric circles in S7 includes the following operations:

[0411] S7.1: For a muscle fiber i that can generate an effective diffusion circle, for the geometric center of the muscle fiber (x i ,y i Generate diffusion circles with three radii: r1 = 1.5·D avg r² = 2.5·D avg r3 = 3.5·D avg , where D avg The global average equivalent diameter of the current type of muscle fiber; or generate three diffusion circles of adjacent layers, representing muscle fibers containing only 1, 1-2, and 1-3 adjacent layers respectively. In this case, r1 = 1, r2 = 2, r3 = 3.

[0412] S7.2: Using the previously calculated scale vectors (taking radian scale as an example, area scale is similar), combine them into a scale matrix:

[0413]

[0414] S7.3: First-order one-step difference calculation: Taking the k-th type of muscle fiber as an example, using three points... Calculate the second-order difference. First, calculate the one-step first-order difference between adjacent points:

[0415]

[0416] S7.4: The second difference is the difference of the first difference:

[0417] S7.5: First-order two-step difference calculation:

[0418] By calculating the arc ratio and area ratio of various types of muscle fibers using diffusion circles of different radii, the composition and distribution patterns of the tissue structure surrounding the muscle fibers, as well as the trends in these changes, can be quantified. (Second-order difference Δ) (2) and first-order difference The combined description of this feature as a function of distance reflects the characteristic pattern of the environment surrounding the muscle fiber and how it changes with distance.

[0419] Here is a specific application example:

[0420] Spatial microenvironment extraction, spatial microenvironment feature calculation, and muscle fiber interaction feature index calculation were performed on SDH enzyme-stained images of skeletal muscle from over 40 sites in Rongchang Large White pigs (abbreviations are shown in Table 4). Table 2 shows the calculation results of some connected component feature indices, and Table 3 shows the calculation results of some feature indices of muscle fiber centrality and diffusion circle. Table 4 shows the English abbreviations for skeletal muscle sites, and Table 5 shows the English abbreviations for muscle fiber feature indices.

[0421] Figure 7 For the clustering demonstration of connected components and diffusion circles, the two images show that each category (A1, A2, A3, B1, B2, B3) reflects significant similarity among connected components, with low similarity between different categories. In category C1, muscle fibers within the diffusion circle are distributed on one side of the circle, while in category C2, muscle fibers are distributed around the center. D indicates the clustering effect of the diffusion circles, showing high similarity (greater than 70%) among the diffusion elements within each category. Therefore, this invention can effectively cluster connected components with similar characteristics and diffusion circles with similar muscle fiber distribution patterns, and distinguish between different categories.

[0422] Figure 8This study reflects the differences in connected component indices among different pig breeds. Taking subgraphs A1, A2, and A3 as examples, the mean (A1) and median (A2) of the connected component rectangularity of Large White pigs are significantly greater than those of Rongchang pigs (effect size +++) (statistically significant***), but the variance of their rectangularity is significantly smaller. This indicates that the shape of the connected components of Large White muscle fibers is closer to a rectangle, with lower shape variability. Therefore, our method can detect differences in indices between pig breeds. These differences in muscle fiber phenotypes may be related to differences in meat quality and mobility among pig breeds.

[0423] Figures 9-10 To differentiate the connected component indices across different body parts, Figure 9 Taking subgraphs (II)A1 and (II)A2 as examples, the mean (A1) and median (A2) rectangularity of the muscle fiber connected regions in the hind limbs of Rongchang pigs are significantly larger than those in the forelimbs, indicating that the shape of the muscle fiber connected regions in the hind limbs of Rongchang pigs is more regular. Figure 9 Taking (I)B2 subgraph as an example, the clustering coefficient of the type I muscle fiber connectivity domain in the forelimbs of the Large White pig is significantly lower than that in the trunk and hindlimbs, reflecting that the connections between muscle fibers in the type I muscle fiber connectivity domains of the forelimbs are more dispersed and sparse. Figure 10 It can be seen that there are significant differences in the ratio of type I muscle fiber contact circumference between Large White pigs and Rongchang pigs in the forelimbs, trunk, and hindlimbs. Specifically, the Large White pigs have a higher degree of contact between type I muscle fibers in the trunk compared to the forelimbs and hindlimbs; while in Rongchang pigs, the contact between type I muscle fibers is higher in the forelimbs. This demonstrates that our method can identify differences in indicators across different body parts of the same pig breed. These differences in muscle fiber phenotypes may be related to differences in meat quality and locomotion between different body parts.

[0424] Figure 11 This is a heatmap showing the percentage difference in the mean of the muscle fiber clustering coefficient and the percentage difference in phenotypic differences between different skeletal muscle sites, calculated using this scheme for sampling different individuals within the same skeletal muscle region. For example, overall, we can see that the differences are smaller in the diagonal region, meaning there is less phenotypic difference within the same skeletal muscle, while the differences are larger in the off-diagonal region, indicating phenotypic variability between different skeletal muscle sites. This demonstrates that the phenotypic index designed in this scheme is conservative within a site and can be used as a representative feature, while exhibiting significant variability between sites and can be used as a differential feature. At the same time, some skeletal muscle sites with extremely large differences can also be observed, such as the deltoid-scapular region (DSP) and the longissimus dorsi (LDM) region, with a difference of 0.62.

[0425] Figure 12 and Figure 13 For data sampled from the same skeletal muscle region (different regions, different individuals) and data sampled from different skeletal muscle regions, their connected component indices ( Figure 12 The contact ratio between the muscle fibers and the muscle fibers ( Figure 13 The mean differences (between and within) of the indicators are the average values ​​of the indicators in the three body regions (divided into forelimbs, hindlimbs, and trunk). It can be seen that, regardless of the forelimbs, hindlimbs, or trunk, the average differences of the characteristic indicators of the connected regions and the muscle fiber contact ratio indicators of different skeletal muscle sites are generally greater than or much greater than the differences of the connected region indicators of the same skeletal muscle sites, further reflecting the representativeness of the characteristic indicators designed in this scheme.

[0426] Depend on Figure 14 It can be seen that the area of ​​type I muscle fibers in various skeletal muscle regions of Rongchang pigs and Large White pigs is significantly positively correlated with their interaction indices (the proportion of type I muscle fiber area within a 2.5-radius diffusion circle (2.5AED-TMCR), the degree of ingress and egress contact with type I muscle fibers, the degree centrality (DegCent) in the connected region, and the betweenness centrality (BetCent)). Overall, this reflects that type I muscle fibers with a higher degree of contact tend to have a larger cross-sectional area. This may mean that there is a regulatory effect promoting muscle fiber hypertrophy between type I muscle fibers in the two pig breeds.

[0427] Depend on Figure 15 and Figure 16 It can be seen that the pH value and conductivity of the longest muscle mass in the back of the Large White pig are significantly correlated with the number of muscle fibers, cluster coefficient, rectangularity, and aspect ratio of the connected domain characteristics of type I muscle fibers. At the same time, they are significantly correlated with the spatial adjacency (in-out contact ratio) between type I muscle fibers and the proportion of type I muscle fibers (variance) in the diffusion circle.

[0428] Depend on Figure 17 and Figure 18 It can be seen that the intramuscular fat content of the longissimus dorsi and psoas major muscles of Large White pigs and Landrace pigs is significantly correlated with the number of muscle fibers, cluster coefficient, rectangularity, aspect ratio, and other indicators of the connected domain characteristics of type I muscle fibers; at the same time, it is significantly correlated with the spatial adjacency degree (in-out contact ratio) between type I muscle fibers and the proportion of type I muscle fibers (variance) in the diffusion circle.

[0429] In summary, it can be seen that the method system provided by this invention can effectively measure the spatial microenvironment characteristics and muscle fiber interaction characteristics of skeletal muscle fibers, and can effectively support downstream differential analysis. It also shows a high correlation with meat quality indicators and is expected to combine molecular-level information to analyze the causes of meat quality.

[0430] The technical features of the above embodiments can be combined in any way. For the sake of brevity, not all possible combinations of the technical features in the above embodiments are described. However, as long as there is no contradiction in the combination of these technical features, they should be considered to be within the scope of this specification.

Claims

1. A method for quantifying the microenvironment of porcine muscle fibers and their interactive features, characterized by, The method comprises the following steps: S1: image information extraction, based on the muscle fiber cross-section imaging figure which has been matched with muscle fiber ID, category, contour coordinate set and geometric center coordinate, the center distance matrix and the shortest distance matrix between muscle fibers are calculated, and the short diameter data of each muscle fiber is calculated by using the minimum circumscribed rectangle to construct an adjacency matrix; S2: diffusion circle microenvironment extraction, starting from the geometric center of each muscle fiber, 1.5, 2.5 and 3.5 times of the average diameter of all muscle fibers are taken as the radius to intercept three diffusion circles; starting from the geometric center of each muscle fiber, the adjacency relationship of muscle fibers is identified based on the distance matrix of muscle fibers, and three adjacent diffusion circles are generated within the range of only the first adjacent circle layer, the first and second adjacent circle layers, and the first, second and third adjacent circle layers; S3: connected domain microenvironment extraction, based on the distance matrix of muscle fibers, the connected domain formed by the adjacency relationship of muscle fibers of the same type is identified; S4: muscle fiber interaction feature measurement, the contact degree feature of muscle fibers and their respective types of muscle fibers is measured, the multi-circle layer adjacency information in the diffusion circle and the center degree in the connected domain are measured; S5: microenvironment feature index measurement, the size index and shape index of the connected domain are measured, the composition ratio index and distribution mode index of each type of muscle fiber in the diffusion circle are measured, and the statistical amount of muscle fiber attributes in the two microenvironments is measured; S6: microenvironment classification, according to the rotation overlap degree and the weighted Euclidean distance between microenvironment feature vectors, weighted k-means and PCA k-means clustering are used for similarity clustering of diffusion circles and connected domains; S7: change trend measurement of diffusion circle features with increasing radius, the proportion vector index in the three adjacent diffusion circles formed by each muscle fiber according to 1.5, 2.5 and 3.5 times of the average diameter of all muscle fibers, and the first adjacent circle layer, the first and second adjacent circle layers, and the first, second and third adjacent circle layers are calculated.

2. The method of claim 1, wherein the microenvironment of porcine muscle fiber composition and its interactive feature quantification is characterized by, The distance matrix construction in S1 comprises the following operations: The center distance between any two muscle fibers in the computed image is calculated, along with the shortest distance, where the center distance between the ith muscle fiber and the jth muscle fiber is denoted by cd ij is denoted by cd Shortest distance with md ij denotes: where Ω i and Ω j represent the set of profile coordinates of the ith muscle fiber and the jth muscle fiber, respectively, the distance matrix CD= {cd ij , the distance matrix MD= {md ij}.

3. The method of claim 1, wherein the microenvironment of porcine muscle fiber composition and its interactive feature quantification is characterized by, The diffusion circle microenvironment extraction in S2 comprises the following operations: S2.1: for all muscle fibers, construct an adjacency matrix, initialize the value of the adjacency matrix ADJ as the distance matrix MD, mark the elements with a value less than 1 / 2 of the minimum length of the short diameter of all muscle fibers as MSD, and set the value to 1 to complete the construction of the adjacency matrix ADJ; The adjacency matrix ADJ has elements satisfying: where md ij is an element of the matrix of shortest distances MD = {md ij} between the i-th and j-th muscle fiber profile, and MSD is the minimum value of the short diameter of all muscle fibers in the picture. S2.2: the first type of diffusion circle is uniquely represented by the microenvironment center point and the microenvironment absolute radius, wherein the microenvironment center point is the geometric center of a muscle fiber, matching the muscle fiber ID, and the circle with a contour exceeding the image boundary is excluded; S2.3: The second type of diffusion circle is uniquely represented by the microenvironment center point and the microenvironment adjacency radius, which refers to the radius defined by the distance of muscle fibers in the adjacency matrix. If muscle fiber i needs to move k steps to reach muscle fiber j through adjacency, their distance is k. For each valid diffusion circle under the three radius lengths of adjacency distances 1, 2, and 3, obtain the muscle fiber ID of the circle center, the circle center coordinates, and the adjacency radius of the diffusion circle as the identifier of this microenvironment. Calculate the contained muscle fiber ID. Take the diffusion circle with an adjacency radius of k for the ith muscle fiber as an example. The contained muscle fiber ID is calculated as follows: find(ADJ(i,:)+ADJ 2 (i,:)+…+ADJ k (i,:)) where ADJ k (i,:) represents the i-th row of the adjacency matrix k times, the find function returns the position of the non-zero elements in it.

4. The method of claim 1, wherein the microenvironment of porcine muscle fiber composition and its interactive feature quantification is characterized by, The S3 microenvironment extraction of connected domains includes the following operations: S3.1: Construct an adjacency matrix ADJ for all muscle fibers of the k-th type (k) ADJ (k) = ADJ(N k ,N k ), where N k is the ID of all muscle fibers of the k-th type in the image; S3.2: Connected domain identification: identify the connected regions for each type of muscle fiber using the breadth-first algorithm.

5. The method of claim 1, wherein the microenvironment of porcine muscle fiber composition and its interactive feature quantification is characterized by, The S4 measurement of muscle fiber interaction features includes the following operations: S4.1: Calculate the contact perimeter ratio of each muscle fiber with each type of muscle fiber; S4.2: Calculate the weighted sum of the contact ratio of each muscle fiber with the surrounding muscle fibers, and select specific weights to reflect the specific level of interaction; S4.3: Calculate the information within the diffusion circle centered on the muscle fiber. Take the kth muscle fiber as an example. The diffusion circle information centered on it is condensed into a matrix consisting of distance, angle, category, and area occupied; S4.4: Measure the interaction intensity of muscle fibers with each type of muscle fiber in the specified local area; S4.5: Measure the position center degree of muscle fibers in the connected domain formed by the same type of muscle fiber.

6. The method of claim 1, wherein the microenvironment of porcine muscle fiber composition and its interactive feature quantification is characterized by, The S5 microenvironment feature index measurement includes the following operations: S5.1: Diffusion circle microenvironment feature quantification, including: S5.1.1: The proportion of the total area of each type of muscle fiber in the diffusion circle. If the category membership is obtained, it is classified according to the 0.5 threshold. Take the i-th muscle fiber as an example, and the diffusion circle with radius or adjacent circle layer r is denoted as ARatio i r ; S5.1.2: Arc proportion of each type of muscle fiber within the diffusion circle; S5.1.3: adjacent information matrix of the ith muscle fiber radius or adjacent shell layer of the preceding configuration r S5.1.4: Interaction attribute statistics of each type of muscle fiber; S5.2: Calculate the autocorrelation Moran coefficient of the index value of the muscle fiber within the diffusion circle; S5.3: Calculate the local distance matrix entropy value and angle matrix entropy value of each type of muscle fiber within the diffusion circle; S5.4: Measure the aggregation degree of each type of muscle fiber within the diffusion circle; S5.5: Connected domain microenvironment feature quantification, including the following steps: S5.5.1: Quantify size indicators: area and number of muscle fibers; S5.5.2: Quantify shape indicators: aspect ratio, circularity, filling degree, rectangularity S5.5.3: Calculate the cluster coefficient.

7. The method of claim 1, wherein the microenvironment of porcine muscle fiber composition and its interactive feature quantification is characterized by, The S6 microenvironment classification includes the following operations: S6.1: For absolute diffusion circles, use the sparse matrix clustering method; S6.2: For connected domains, use weighted k-means and PCA-Kmeans clustering; S6.3: For the valid adjacent diffusion circle, the selection of the features is the proportion of the total area of each type of muscle fiber, the proportion of the arc of each type of cell, the local distance matrix entropy value and the angle matrix entropy value of each type of muscle fiber, and the aggregation degree of each type of muscle fiber, which are combined into a feature vector instead of a matrix. If the adjacent radius exceeds 1, for example, k, the adjacent diffusion circle features are calculated for the adjacent radii i = 1, 2, …, k, respectively, and then these features are spliced into a feature vector. Next, according to the feature vectors, the PCA-Kmeans is used for clustering.

8. The method of claim 1, wherein the microenvironment of porcine muscle fiber composition and its interactive feature quantification is characterized by, The S7 change trend measurement of the diffusion circle features with the increase of the radius or the adjacent circle layer includes the following operations: S7.1 : For a certain muscle fiber i for which an effective diffusion circle can be generated, generate diffusion circles of three radii for the geometric center (x i ,y i ) of the muscle fiber: r1 = 1.5 · D avg , r2 = 2.5 · D avg , r3 = 3.5 · D avg , where D avg is the global average equivalent diameter of the muscle fiber of the current class. S7.2: Using the previously calculated proportion vector, a proportion matrix is combined: S7.3: First-order one-step difference calculation: Taking the k-th type of muscle fiber as an example, a three-point To calculate the second-order difference, first calculate the one-step first-order difference between adjacent points: S7.4: Difference of second order difference to first order difference: S7.5: First order two-step difference calculation:

9. The method according to claim 7, wherein the microenvironment of porcine muscle fiber constitution and its interactive characteristic quantification is characterized in that, The weighted k-means includes the following steps: (1) Data reading and feature normalization: The pandas library is used to read the data in the connected domain information document, and the features that need to be clustered are features = ['circularity', 'aspect ratio', 'filling degree','rectangularity', 'total area', 'clustering coefficient'], which are normalized using MinMaxScaler respectively. The normalization formula is: where x is the original feature value, i.e. the specific value of a sample in the data set on a certain feature; (2) Weight setting: custom weight list w * = [0.1, 0.1, 0.1, 0.1, 0.5, 0.1], because the total area is more intuitive, so the weight setting is the largest; (3) The elbow rule is used to determine the optimal number of clusters optimal_clusters, and the maximum number of clusters num_clusters is set to 20; For each k cluster , solve the objective function: where, denotes different clustering partitioning schemes, the goal is to find the optimal partitioning scheme that minimizes the objective function value, double summation: outermost sum over all clusters, middle layer sum over all samples within each cluster and the squared distance to the cluster center, weighted Euclidean distance: wherein, is the normalized feature vector of the jth connected component, u i is the center (mean) of the connected component feature vectors in the ith cluster, is the weight of the kth feature. The second derivative is calculated: SSE''(k) ≈ SSE(k+1)-2SSE(k)+SSE(k-1), where k ∈ [2, num_clusters-1], k = 1 cannot calculate k-1 = 0; k = num_clusters, cannot calculate k+1 = num_clusters+1 (out of the initial set range); The optimal number of clusters is:

10. The method of claim 7, wherein the microenvironment of porcine muscle fiber composition and its interactive characteristic quantification is characterized by, The PCA-Kmeans clustering includes the following steps: (1) Data reading and feature screening: Taking connected domain as an example, the connected domain data is read, 6 original features features = ['circularity', 'aspect ratio', 'filling degree','rectangularity', 'total area', 'clustering coefficient'] are selected for analysis, an original feature matrix X is constructed, the number of connected domains of an image is n, the number of features p = 6, and the original feature matrix X = [x ij ]∈R n×p , wherein x ij represents the jth feature value of the ith connected domain; (2) Z-score standardization: The original features are standardized; wherein (3) Principal component analysis: First, the eigenvalues and eigenvectors of the covariance matrix of the standardized data are calculated; Covariance matrix: Eigenvalue and eigenvector computation: solving ∑v = λv k = λ k v k where λ k is the eigenvalue and v k is the corresponding unit eigenvector; Then, the explained variance ratio is calculated, and the explained variance ratio of a single principal component is: wherein λ k is the eigenvalue corresponding to the kth principal component (kth eigenvalue), and p is the total number of original features, the sum of all principal component eigenvalues; Cumulative explained variance proportion: where k represents the first k principal components (k ≤ p); The cumulative explained variance ratio represents the proportion of the total variance of the original data explained by the first k principal components; The number of principal components is determined, and when n ≥ p, k is selected to make the cumulative variance ratio ≥ 85%; when n < p, k = min(n-1, p) and the cumulative variance ratio ≥ 80%; The principal component scores are calculated, the standardized data is projected into the principal component space, and a principal component score matrix Z after dimension reduction is obtained: Z = X scaled • V, where V = [v 1, v 2, v 3, …,v k ]∈R p×k ; (4) Directly use the principal component scores in Z to cluster: The elbow rule is used to determine the optimal number of clusters optimal_clusters, and the maximum number of clusters num_clusters is set to 20; For each k cluster , solve the objective function: to determine a clustering solution, where C i is the set of indices contained in the i-th cluster, u i is its center vector in the Z matrix, i.e. the mean of the corresponding principal component score vectors; The second derivative is calculated: SSE''(k) ≈ SSE(k+1)-2SSE(k)+SSE(k-1), Wherein, k ∈ [2, num_clusters-1], k = 1 cannot calculate k-1 = 0; k = num_clusters, cannot calculate k+1 = num_clusters+1 (out of the initial set range); Optimal cluster number:

Images