Literature DB >> 30181964

An Optimized Framework for Cancer Prediction Using Immunosignature.

Fatemeh Safaei Firouzabadi1, Alireza Vard2, Mohammadreza Sehhati2, Mohammadreza Mohebian3.   

Abstract

BACKGROUND: Cancer is a complex disease which can engages the immune system of the patient. In this regard, determination of distinct immunosignatures for various cancers has received increasing interest recently. However, prediction accuracy and reproducibility of the computational methods are limited. In this article, we introduce a robust method for predicting eight types of cancers including astrocytoma, breast cancer, multiple myeloma, lung cancer, oligodendroglia, ovarian cancer, advanced pancreatic cancer, and Ewing sarcoma.
METHODS: In the proposed scheme, at first, the database is normalized with a dictionary of normalization methods that are combined with particle swarm optimization (PSO) for selecting the best normalization method for each feature. Then, statistical feature selection methods are used to separate discriminative features and they were further improved by PSO with appropriate weights as the inputs of the classification system. Finally, the support vector machines, decision tree, and multilayer perceptron neural network were used as classifiers.
RESULTS: The performance of the hybrid predictor was assessed using the holdout method. According to this method, the minimum sensitivity, specificity, precision, and accuracy of the proposed algorithm were 92.4 ± 1.1, 99.1 ± 1.1, 90.6 ± 2.1, and 98.3 ± 1.0, respectively, among the three types of classification that are used in our algorithm.
CONCLUSION: The proposed algorithm considers all the circumstances and works with each feature in its special way. Thus, the proposed algorithm can be used as a promising framework for cancer prediction with immunosignature.

Entities:  

Keywords:  Cancer; feature selection; immunosignature; normalization

Year:  2018        PMID: 30181964      PMCID: PMC6116316          DOI: 10.4103/jmss.JMSS_2_18

Source DB:  PubMed          Journal:  J Med Signals Sens        ISSN: 2228-7477


Introduction

When antibodies circulate in the blood, they can connect to a large microarray of randomized sequence peptides.[1] An “immunosignature” is a pattern of random sequence peptides, which is obtained with a blood sample test.[2] Neoantigens are produced by cancer release native proteins and biomolecules that are not encountered by the immune system. Therefore, the change in the regulation of the gene expression and proteins in cells can be considered as a sign of cancer.[3] However, there is a slight overlap of the signatures among the cancers that resulted in a loss of specificity in distinguishing between cancers using immunosignatures. To resolve this, peptides were determined to be statistically significant in the cancer signatures using more stringent selection processes. In the recent years, various methods have been introduced in the literature for predicting cancer with peptides and proteomic datasets. Zhang et al.[4] for classifying ten types of cancers used the protein expression profiles. They used the minimum redundancy maximum relevance and the incremental feature selection methods in order to select 23 out of 187 proteins as the inputs of the sequential minimal optimization classifier. Using 23 proteins, they classified with Matthews Correlation Coefficient (MCC) of 0.936 on an independent test set. Kaddi and Wang[5] used the proteomic and transcriptomic data to predict the early stage of cancer in head-and-neck squamous cell carcinoma. They proposed a filter and wrapper method for feature selection and employed the individual binary classification accompaniment with the ensemble classification methods. Stafford et al.[2] used ANOVA and t-test for approximately 10,000 peptide sequences. They proposed a novel feature selection method and the naive Bayes, linear discriminant analysis, and support vector machine (SVM) for classification in two libraries. An average accuracy of 98% and an average sensitivity of 89% were reported.[6] Nguyen et al.[7] used five filter-based feature selection methods including t-test, Wilcoxon, entropy, signal-to-noise ratio, and receiver operating characteristic curve. Moreover, they employed an analytic hierarchy process which is a multi-criteria decision analysis based on the type-2 fuzzy method for the classification of cancer microarray gene expression profiles. Based on the proposed classification method, the achieved accuracy was 95.24%, considering t-test as the feature selection method. In the proposed framework, we combined the meta-heuristic population-based optimization with feature selection and normalization methods to improve the accuracy and efficiency of the classification algorithms for classifying 12 different cancer types. Briefly, the particle swarm optimization (PSO) method at first filters the significant peptides. Next, it selects the best method of normalization from the dictionary for one feature and chooses weights for features that are selected with the statistical feature selection process. Then, the selected features apply to classification. In our study, three types of classification, including SVM, multilayer perceptron (MLP), and decision tree (DT), are used and compared with each other. The rest of the article is organized as follows: in the next section, information about the datasets and proposed methods which are used in this study are presented. The results of the proposed method and the discussion about them are provided in the Results and Discussion Sections and finally, the article is concluded in the Conclusion Section.

Materials and Methods

The dataset

In this study, we used a public immunosignature peptide microarray dataset (arrays of 10,371 peptides) which consists of 1516 patients for 12 different cancer types, 2 infectious diseases, and healthy controls.[7] The Gene Expression Omnibus series code of this dataset is GSE52581, which is publicly available, and other types of this kind of dataset are used in recent studies carried out on cancer.[289] The features of this dataset are not normalized and the dataset includes cancers such as astrocytoma, breast cancer, multiple myeloma, lung cancer, oligodendroglia, ovarian cancer, advanced pancreatic cancer, and Ewing sarcoma. We removed the data related to the infectious diseases and used 1292 subjects. Thus, the dataset consists of 1292 columns and 10,371 rows and, in this article, rows refer to mean peptides and columns refer to samples. More details about the used dataset are described in Table 1.
Table 1

Basic information of patients per cancer in the utilized dataset

Basic information of patients per cancer in the utilized dataset

The proposed algorithm

The structure of the proposed algorithm is depicted in Figure 1. As shown in this figure, at first, features are normalized with different methods that are selected by PSO from a dictionary of methods. Then, a statistical feature selection method is used to identify significant discriminative features. The selected features are assigned weights by PSO. The goodness of fit of PSO can be measured by F1-score of classifier on the training set (Eq. 1).
Figure 1

The structure of the proposed algorithm

The structure of the proposed algorithm Goodness of fit = FScoretrain (1) After selecting suitable features, classification methods are used to classify the types of cancers. In this study, we utilized three multi-class classification methods including multi-class SVM, DT, and multilayer perceptron. The weights of features are estimated using PSO during learning of classifiers. The algorithm stops if no remarkable improvement is seen in the objective function or the maximum number of iterations (set to 50 in our study) is reached. In the following section, the methods and algorithms, employed in different parts of the proposed scheme, are described in detail.

Normalization

According to Liu et al.,[10] a major bioinformatics challenge is the normalization of the data. Normalization is a method that puts data in a similar domain when they are not in one domain. In other words, a data miner may encounter situations in which features in data include quantities in different domains. These features with large quantities may have higher impacts on cost function compared to features with low quantities. This issue is resolved by normalization of features such that their values are put into one domain.[11] Thus, if each feature normalizes properly, classification would be applied more effectively in feature space compared to using classifying without normalization.[12] In addition, it changes the characteristics of the underlying probability distributions.[13] A careful analysis of the geometry of feature space suggests a modification on normalization procedure that works on each feature separately.[14] Therefore, in this study, a dictionary of normalization methods is used which gives the algorithm the ability to select the best normalization method for each feature irrelevant from other features. These normalization methods are selected from the previous studies used on peptides, proteomics, and other microarray datasets.[12151617] These methods are described in the following subsection. Global median centering In this method, the median value of peptides is simply subtracted from each peptide. While this method needs only one sample and does not require other data column, it has bias when the number of peptides is low (<100).[12] Tukey (median polish) This method is a nonparametric analysis of variance.[12] The column medians and row medians are subtracted until a defined criterion is satisfied. This criterion could be the number of iteration or a specific value of the proportional reduction in the sum of absolute residual.[18] If the number of samples and peptides is reasonably large, the approach will be robust.[11] In this study, the termination criterion of the algorithm is empirically defined as 400 iterations. Loading control method This method uses the global median-centering approach at the first step for each peptide. Then, global median-centering approach is used for each sample rather than each peptide.[12] Robust Z-score This method is robust like standard score and is formulated with Eq. 2.[19] Where X’ij represents the ith peptide and the jth sample. In addition, MADi and Mediani denote median and median absolute deviation, respectively. According to the formulation, this method is not sensitive to outliers and therefore is appropriate for microarray data.[20] Invariable method This method is well known for microarray data analysis in the literature[122122] since it can eliminate systematic variation in microarray data. In this method, peptides are ranked, and then a peptide with the highest rank is discarded. This process is repeated until the remaining number of peptides reach a predetermined value (1000 in this article). Then, 25% of the highest and the lowest ranked data are removed and the average of each peptide creates a virtual reference. At the end, each sample is normalized to the virtual reference by the MA-plot approach.[12] In the MA-plot approach, the difference between each sample and the reference sample in the logarithm of base 2 is plotted against the mean of each sample and the reference sample. Then, the normalized values are generated using residuals of fit consequently. Modified Z-score method In this method, log2 is performed on the whole data, firstly. Then, the standard score normalization is applied to each sample and each peptide subsequently. In the final step, the arctangent function is applied to the data.[2324] This method has had a great performance on gene expression data comparatively.[23242526]

Statistical feature selection

Feature selection approaches can be classified into three categories including wrapper, filter, and embedded methods.[27] The filter method is used in the preprocessing section and works independently from classifier.[28] On the other hand, other two methods are used during classification. The wrapper method evaluates the combination of features by formulating a problem and searches the problem space for the best features.[29] This method tests the entire feature subsets.[30] Finally, the embedded method evaluates the accuracy of the classifier for predicting the best features with searching that is guided by a learning classification process. This characteristic of the embedded feature selection method makes it robust against overfitting.[27] According to the characteristics of the peptides’ datasets that have only interval data types, the Kolmogorov–Smirnov (KS) test was appropriate.[31] The KS evaluates the maximum absolute difference in the overall distribution of the two groups (cancer or noncancer). Then, if a feature was significant, the independent-sample t-test was used to identify statistically discriminative normally distributed features.[32] Otherwise, the Mann–Whitney test is performed to check whether two independent samples are significantly different or not.[33] Finally, if numbers of selected features are >50, the algorithm selects the 50 highest ranks in t-test and then in Mann–Whitney test. Thus, the 50 discriminative features are selected consequently and prepared as the inputs of the classification system which is combined with PSO for weighing them as discussed in the following subsection.

Particle swarm optimization

PSO is a meta-heuristic stochastic evolutionary population-based optimization algorithm that is inspired by birds.[34] PSO can solve a range of difficult optimization problems, but it has shown a faster convergence rate compared to other evolutionary algorithms.[3536] This algorithm is used in many different computational biology fields such as modeling in biology,[37] feature selection in gene expression data,[3839] DNA sequence encoding,[40] and breast cancer recurrence prediction.[6] The PSO algorithm starts with random solutions that are called particle positions. Each particle has a velocity and a position that help it to search the whole problem space. The position of a particle is updated according to three parameters: its previous speed, the best position visited by particle, and the best position of the neighborhood (Eq. 4). where n is the iteration number; C1 and C2 are the learning factor coefficients that usually set to 2; i is the particle number, r1, and r2 are random numbers that are uniformly distributed in (0, 1); w is the inertia weight, where a large number of it shows exploration while a small number of it denotes exploitation[41] Thus, the inertia coefficient was set to 1.00 at the first iteration and was linearly decreased with the damping coefficient of 0.99 at each iteration. Furthermore, the objective function is maximizing the goodness of fit that is mentioned in Eq. 1. In our study, PSO algorithm has three main tasks. First, it selects probable peptides. For this purpose, 400 particles are considered in which each particle represents the index of one peptide and it can be an integer number between 1 and 10,371. Second, the PSO algorithm chooses the normalization method for each feature that was selected in the previous step. In these steps, the rounded value of the particle position is selected as an index of the normalization method. It is necessary to say, only one peptide which is normalized individually or with accompaniment of other peptides goes to the next step; for example, if the loading control method is selected for one peptide, the peptide will be normalized by this method with accompaniment of other peptides (because other peptides are used in normalization formulation). However, other peptides will be normalized with their own selected method and then move to the next step. Third, the selected features of the statistical feature selection method are weighed by the PSO algorithm that weighs a real number between 0 and 1. In a nutshell, 400 particles for initial filtering of features, 400 particles in selecting normalization method, and 50 particles for weighing features were used. In total, 850 variables should be considered by PSO, which is an appropriate algorithm for solving high-dimensional problems.[4243] After features are selected by statistical feature selection and weighed by PSO, they are used as inputs of the classifiers.

Classification

In the proposed method, three classification methods such as MLP, DT, and SVM are used. These methods were used in previous similar studies.[6444546] The SVM considers a set of the hyperplane in a high-dimensional space.[47] This algorithm has been widely utilized in peptides’ datasets and other relevant fields.[64849] The radial basis function (RBF) kernels were used in this study and it is paramount importance that the soft-margin parameter and the radius of the RBF kernel should be set appropriately since they do not cause poor classification results. We used the method proposed by Wu and Wang[50] for this purpose. MLP is an improved version of the standard linear perceptron method and can be used for classification of nonlinearly separable data.[51] It is a famous machine learning approach in a variety of the computational biology field such as prediction of protein stability, prognosis DNA methylation biomarkers in ovarian cancer, and encoding aminoacids.[525354] In this study, the MLP with one hidden layer, 20 neurons, and the sigmoid activation function are used because they seem suitable for prediction of cancer, according to the literature.[555657] The number of neurons is calculated empirically until increasing the number of neurons up has no effect on performance. The DT classifier creates a structure of tree for modeling. We used C4.5 DT classifier in the proposed model, which is an entropy-based algorithm that can handles continuous attributes.[5859] Furthermore, this method is widely proposed in predicting cancer, predicting specific target of peptides, and analyzing microarray dataset.[60616263]

Validation procedure

The performance of the proposed algorithm was evaluated by the hold-out method in terms of sensitivity, specificity, precision, accuracy, F1-score, and MCC, which are defined in Table 2.
Table 2

Validation parameters

Validation parameters The sensitivity and specificity are highly dependent on the prevalence of the diseases. On the other hand, a reliable diagnostic system has sensitivity and specificity more than 80% (the minimum statistical power of 80%) and 95% (the maximum Type I error of 0.05), respectively.[6465] Thus, a conservative method should satisfy both parameters.

Results and Discussion

Leave-one-out cross-validation often underestimates error and leads to overfitting.[2] Therefore, in this study, we used hold-out validation and the results of the proposed method on the test set are shown in Tables 3–5. In this study, the hold-out method is used with 70% data as a training set and 30% data as a test set.[6] The 70% and 30% of data were chosen randomly for training and test sets in each time of running the whole algorithm. The complete procedure of data processing has been done 7 times and the results that are shown in Tables 3–5 are mean of 7 times of running. These results indicate that the procedure with the proposed methods is robust. The final classification system is a multiclass system in which its parameters are calculated based on systematic analysis of multiclass classification approach[66] which are shown in Tables 3–5.
Table 3

The average of five holdout performance estimates for multilayer perceptron method

Table 5

The average of five holdout performance estimates for decision tree method

The average of five holdout performance estimates for multilayer perceptron method The average of five holdout performance estimates for support vector machine method The average of five holdout performance estimates for decision tree method The F-score on the training set and the test set against PSO iteration is depicted in Figure 2.
Figure 2

The value of the fitness function. F-score on the training set is the solid line and the F-score on the test set is the dash-dot line during optimization procedure. The termination criterion was only the maximum number of iterations (i.e., 20) in this plot

The value of the fitness function. F-score on the training set is the solid line and the F-score on the test set is the dash-dot line during optimization procedure. The termination criterion was only the maximum number of iterations (i.e., 20) in this plot The average time measured for validating 454 patients were 0.33 ± 0.07, 0.30 ± 0.05, and 0.34 ± 0.08 second, respectively, for SVM, DT, and MLP. All results were obtained on a computer of Intel Core-i7, 2 GHz CPU with 8 GB of RAM. Due to data imbalance, it seems that overall accuracy is not a suitable fitness measure.[67] As a matter of fact, selecting inappropriate objective function instead of Eq. 1 creates bias toward majority prevalence. The average sensitivity and specificity of the algorithm based on SVM classifier are estimated as 99.8 ± 0.7 and 99.9 ± 0.6 [Table 4], respectively, as the best classifier in this approach. Furthermore, the maximum error type I (α) and II (β) are 0.009 and 0.076, subsequently. Thus, in the proposed method, the average of type I and II errors showed consistency in the results of the algorithm and this method has the capability to be used in clinical applications.
Table 4

The average of five holdout performance estimates for support vector machine method

The proposed method could be called a general framework since it uses a dictionary of normalization's method and optimized feature selection. As a comparison with others’ works which just use t-test or Wilcoxon feature selection with a fixed type normalization method, this method could search more space of the solutions and present a comprehensive answer. In other words, methods which only use one normalization method and feature selection procedure are one of the solutions in the search space of the proposed framework. As the purpose of illustration, Stafford et al.[2] worked on the same dataset and used t-test method as feature selection and global median centering as normalization that gained an average accuracy of 98% and an average sensitivity of 89%, while in the proposed method, the average accuracy was 99.16% and the average sensitivity was 95.87% which revealed the advantages of the proposed method.

Conclusion

This article provides a novel method for predicting cancer with immunosignature. In the proposed method, the PSO algorithm was first used to filter some features. The selected features were refined by the statistical feature selection methods and estimated weights by PSO. The overall feature selection process is performed as part of a learning procedure. Instead of PSO, other meta-heuristic population-based stochastic optimization methods that could deal with discrete (feature and normalization selection) and continuous (feature weights) problems could be used. The performance of the algorithm was dependent on the PSO initialization and classifier tuning. As an instance, choosing proper numbers of hidden layers and neurons in each layer in MLP, kernel's parameters in SVM, and search space for initialing weight on feature selection by PSO procedure could have the paramount effect on results. In the proposed method, the normalization dictionary is independently used for each feature because each normalization method was tested on the statistical selected features; however, we cannot confidently select an appropriate normalization method to reach high performance; therefore, an optimization framework was designed. In a nutshell, the modified Z-score normalization method is selected more than other methods by PSO optimization. To shed light on it, it seems that modified Z-score normalization could map features to the suitable new space, in which discrimination between classes would be made more effective than other subspaces. In the proposed study, different algorithms were analyzed on immunosignature data to give a clear insight into the classification and identification of biological markers for the diagnosis of diseases. It can help to adopt useful approaches to early diagnoses and treatments. More specifically, this study proposes a comprehensive algorithm by presenting different methods of normalization and feature selection using PSO which can help to attain optimal results. The proposed algorithm is promising and can be utilized as a new offline tool in clinical applications. The developed program is available to interested readers upon request. One of the limitations of the current study is that results might have been biased because scant of enough available samples. The sample size must be increased to improve the statistical power in our diagnosis system.[68] Another limitation of the proposed algorithm is that output of the classification system is not fuzzy.[69] It will be useful to report the risk of having a cancer type, which is the focus of our future activity.

Financial support and sponsorship

Nil.

Conflicts of interest

There are no conflicts of interest.

BIOGRAPHIES

Fatemeh Safaei Firouzabadi was born in Tehran, Iran, in 1989. She received her undergraduate degree in Physics Engineering from Islamic Azad University, Science and Research Branch, Tehran, Iran in 2013. In 2018, she succeeded to complete a master program in Biomedical Engineering from Isfahan University of Medical Sciences, Isfahan, Iran. She focused her master thesis on “Finding an appropriate method for extraction of protein biomarkers to classify different cancers using protein microarray”. Her research interests are Bioinformatics, Computational Biology, and Biological Modelling. Email: fatemeh.safaei.f@gmail.com Alireza Vard received his B.Sc. degree in software engineering from University of Isfahan, Isfahan, Iran in 2004, and he got his M.Sc. and Ph.D. degrees in computer engineering from University of Isfahan, in 2007 and 2012 respectively. He is currently an assistant professor at Department of Biomedical Engineering, Isfahan University of Medical Sciences, Isfahan, Iran. His research interests include medical image processing, pattern recognition, design and development of medical data processing software. Email: vard@amt.mui.ac.ir Mohammadreza Sehhati is an assistant professor at the Isfahan University of Medical Sciences, in the Biomedical Engineering Department. He received the BS and MS degrees in biomedical engineering from Shahed University and University of Tehran, Tehran, Iran, respectively. He obtained his PhD at the Isfahan University of Medical Sciences in 2015. His research interests include bioinformatics, machine learning, data mining, and image processing. Email: mr.sehhati@amt.mui.ac.ir Mohammadreza Mohebian studied biomedical engineering at the University of Isfahan. In his B.Sc. project, he looked at “A Hybrid Computer-aided-diagnosis System for Prediction of Breast Cancer Recurrence (HPBCR) Using Optimized Ensemble Learning”. He is graduated with the master degree in biomedical engineering from the University of Isfahan too. In his M.Sc. he worked on “Non-invasive Decoding of the Motoneurons: A Guided Source Separation method based on Convolution Kernel Compensation with Clustered Initial Points”. He ranked in the top 10 percent of the B.Sc. students with the same entrance date and was the first rank in M.Sc. degree in university. Moreover, he developed the first website for using machine learning and data mining in clinical data in Iran. Presently, he is studying on biomedical engineering at electrical and computer department of Saskatchewan University. Email: van_mohebian@yahoo.com
  39 in total

1.  CLIFF: clustering of high-dimensional microarray data via iterative feature filtering using normalized cuts.

Authors:  E P Xing; R M Karp
Journal:  Bioinformatics       Date:  2001       Impact factor: 6.937

Review 2.  A review of feature selection techniques in bioinformatics.

Authors:  Yvan Saeys; Iñaki Inza; Pedro Larrañaga
Journal:  Bioinformatics       Date:  2007-08-24       Impact factor: 6.937

3.  High-throughput protein expression analysis using tissue microarray technology of a large well-characterised series identifies biologically distinct classes of breast cancer confirming recent cDNA expression analyses.

Authors:  Dalia M Abd El-Rehim; Graham Ball; Sarah E Pinder; Emad Rakha; Claire Paish; John F R Robertson; Douglas Macmillan; Roger W Blamey; Ian O Ellis
Journal:  Int J Cancer       Date:  2005-09-01       Impact factor: 7.396

4.  Immunosignature system for diagnosis of cancer.

Authors:  Phillip Stafford; Zbigniew Cichacz; Neal W Woodbury; Stephen Albert Johnston
Journal:  Proc Natl Acad Sci U S A       Date:  2014-07-14       Impact factor: 11.205

5.  Improved normalization of systematic biases affecting ion current measurements in label-free proteomics data.

Authors:  Paul A Rudnick; Xia Wang; Xinjian Yan; Nell Sedransk; Stephen E Stein
Journal:  Mol Cell Proteomics       Date:  2014-02-21       Impact factor: 5.911

6.  Algorithm-driven artifacts in median polish summarization of microarray data.

Authors:  Federico M Giorgi; Anthony M Bolger; Marc Lohse; Bjoern Usadel
Journal:  BMC Bioinformatics       Date:  2010-11-11       Impact factor: 3.169

7.  A Decision Tree Based Classifier to Analyze Human Ovarian Cancer cDNA Microarray Datasets.

Authors:  Meng-Hsiun Tsai; Hsin-Chieh Wang; Guan-Wei Lee; Yi-Chen Lin; Sheng-Hsiung Chiu
Journal:  J Med Syst       Date:  2015-11-03       Impact factor: 4.460

8.  A fuzzy logic diagnosis system for classification of pharyngeal dysphagia.

Authors:  S Suryanarayanan; N P Reddy; E P Canilang
Journal:  Int J Biomed Comput       Date:  1995-03

9.  Enrichment of CD44 in basal-type breast cancer correlates with EMT, cancer stem cell gene profile, and prognosis.

Authors:  Hanxiao Xu; Yijun Tian; Xun Yuan; Yu Liu; Hua Wu; Qian Liu; Gen Sheng Wu; Kongming Wu
Journal:  Onco Targets Ther       Date:  2016-01-21       Impact factor: 4.147

10.  Diagnosis of Ovarian Cancer Using Decision Tree Classification of Mass Spectral Data.

Authors:  Antonia Vlahou; John O. Schorge; Betsy W. Gregory; Robert L. Coleman
Journal:  J Biomed Biotechnol       Date:  2003
View more
  2 in total

1.  Study the Effect of the Risk Factors in the Estimation of the Breast Cancer Risk Score Using Machine Learning.

Authors:  Sam Khozama; Ali Mahmoud Mayya
Journal:  Asian Pac J Cancer Prev       Date:  2021-11-01

2.  A Cost-Effective and Portable Optical Sensor System to Estimate Leaf Nitrogen and Water Contents in Crops.

Authors:  Mohammad Habibullah; Mohammad Reza Mohebian; Raju Soolanayakanahally; Khan A Wahid; Anh Dinh
Journal:  Sensors (Basel)       Date:  2020-03-06       Impact factor: 3.576

  2 in total

北京卡尤迪生物科技股份有限公司 © 2022-2023.