{'fraction_non_alphanumeric': 0.02730148687124328, 'fraction_numerical': 0.031667193925972795, 'mean_word_length': 4.280821917808219, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 37, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Image Bio-markers and Gene Expression Data Correlation Framework for Lung Cancer Radio-genomics Analysis Based on Deep Learning

Background: Precision medicine, a popular treatment strategy, has become increasingly important to the development of targeted therapy. To correlate medical imaging with prognostic and genomic data, researches in radiomics and radiogenomics have provide many pre-deﬁned image features to describe image information quantitatively or qualitatively. However, in previous researches, there are only statistical results which proves high correlation among multi-source medical data, but those can’t give intuitive and visual result. Results: In this paper, a deep learning based radio-genomics framework is provided to construct the linkage from lung tumor images to genomics data and implement generation process in turn, which form a bi-direction framework to map multi-source medical data. The imaging features are extracted from auto-encoder under the condition of genomics data. It can obtain much more relevant features than traditional radio-genomics methods. Finally, we use generative adversarial network to transform genomics data onto tumor images, which gives a cogent result to explain the linkage between them. Conclusions: Our proposed framework provides a deep learning method to do radio-genomics researches more functionally and intuitively.


Background
Currently, the number of lung cancer patients and cancer-related deaths has occupied a vital position of cancer-related deaths worldwide. Over 70% of lung cancer patients can only be diagnosed after the onset of symptoms from advanced local or metastatic disease. Unfortunately, the survival rate is only 50% on the condition that diagnose could be located. It is worst that only less than 20% of these patients are diagnosed with very early stages. In this context, precision medicine is gaining popularity for providing customized or personalized healthcare, and quantitative imaging has been contributing to significant improvement of the diagnosis procedures [1].
In another way, personalized medicine aims to tailor medical care to the individual from view of molecular. High-throughput molecular biology technologies claim to produce biomarkers for disease diagnosis and its prognosis prediction [2,3,4]. However, these biospies generated from heterogeneous lesions cannot completely represent the anatomic and physiologic properties, such as the tumor size, anatomic location, and morphology. In another way, image features extracted from these lesions construct a highly informative pathway to disease diagnosis, treatment planning and clinical analysis. Nevertheless, only few studies have constructed radiogenomics frameworks of these information correlations that integrate the genomic and image data [5,6,7,8].
Traditionally, radiology and image-guided interventional therapy have been used to deal with diagnosis and provide anatomical information. However, all these methods must make the patient's physical wounds painful and take longer to heal. To overcome these shortcomings, radiomics extract image features and sub-visual features from radiological images and use state-of-the-art machine learning techniques to provide unique potential for faster and more accurate lung cancer screening. Therefore, by improving the practice of qualitative and quantitative analysis, it is expected to improve the prognosis prediction and response to some treatments. In addition, some radiogenomics researchers have shown that some image features are even associated with genomic changes in tumor DNA [1]. These characteristics can identify specific changes in biological pathways that in turn affect patient management and health outcomes [11,10].
We develop a radiogenomic framework centered on deep learning (DL) to map image features and genomic data, based on our previous work [12] on the correlation between genomics and images. Conditional autoencoder replace the original model while correlating and more genomic analysis are conducted in this paper. Moreover, inspired by some latent space projection model, we fulfill the transition between tumor image and genomic data.

Radiomics and Radiogenomics Approaches
Traditionally, radiomics and radiogenomics contain four steps: image acquisition, lesion segmentation, feature extraction and model validation. In relevant researches of radiomics and radiogenomics, features are extracted from medical images qualitatively and quantitatively including semantic (prognostic) and numerical (statistic) features. Then statistical methods will give a correlation between them and genomic data using gene-set enrichment analysis (GSEA). These works improve diagnostic and prognostic performance in various oncologic applications and finally promote the development of precision medicine. Coroller et al. build a correlation between image features and clinical data to predict distant metastasis in lung adenocarcinoma [13]. Abdollahi et al. use statistic features to predict sensorineural hearing loss and get a high accuracy [27]. The correlation between genomic data and images can indicate the linkage of gene change and tumor variation. Aerts et al. provide a quantitative method to correlate gene expression profile data and the low level image features, which can give a decision-support in cancer treatment at low cost [14]. Additionally, Gevaert et al. propose a protocol which maps semantic features and genome data, and get a model with an area under the receiver operating characteristic curve (AUC) of 65% or greater [29]. Besides these, many researches have proved that there is linkage between multi-source medical data. However, the pre-defined features used in these researches are not rich and effective enough, such as the sum of pixels with different percentages. This gives us the thought that mapping these data onto DL methods, which will produce multi-level and rich image features for correlation.

Deep Learning for Radiological Analysis
Radiomics and DL are focus points in medical imaging field [13]. Radiomics apply image features for prognosis prediction, which is vital due to its clinical significance. DL has been used in medical imaging analysis tasks, such as in CT, MRI, and PET, because its high precision. It can provide informative information about diagnosis, prognostic data, tumor phenotypes, and the gene-protein signatures in lung cancer treatments and prediction [13]. Some recent progresses are as followed.
As a basic model of Convolutional neural networks (CNNs), autoencoder is a popular encoder model. In medical imaging segmentation task, it can give an anatomical prior which eliminates the burden of providing paired example segmentations [21,19]. As an end-to-end method, autoencoder has been applied for lesion detection and segmentation, pixel repair and prognosis validation [21,22,39], which show the features extracted from autoencoder contain mostly image information.
With the similar structure as autoencoder, the U-Net architecture achieves excellent performance on different biomedical applications, including tumor segmentation and CT reconstruction [17,18]. As most popular method of segmentation field, U-net shows its robustness and effectiveness, so we will use it to obtain the lesion region.
Recently, There are many applications in several medicals imaging fields using GAN [20] frameworks, including creating translation of label-to-image, mask-toimage or medical cross modality [24,25,26]. These generation models can output the images we need according to the input data. However, Conditional Generative Adversarial Nets (CGAN) utilizes extra inputs for fusing more information [23].
It produces more precise result compared to GAN due to the extra condition y in generation process: It has been applied in many medical fields, including organ segmentation [34], lesion generation [33], etc.. These researchers indicate CGAN architecture can generate expected medical images.
Considering the distribution of images of different categories, Bao J. et al. propose a framework CVAE-GAN [35] based on VAE [36] different from GAN. This model projects images of different categories to different latent space which makes generation of images easy while giving a specific category. Inspired by this work, we will propose similar method to visualize tumors using genomics data.

Methods
In this part, we will introduce our methods of the whole framework in Fig. 1 for details. First, an U-net based segmentation method is applied to get TR from original CT image. Then, we use autoencoder to encode the images under the condition of gene. The image features are extracted from different levels of encoder. A series of analysis experiments, including prognosis and GSEA, are applied to these features, prognostic data and genes to prove the correlation among these multi-source data. Finally, a modified CVAE-GAN transforms gene to corresponding TR and give an intuitive result.

Tumor Detection and Segmentation
In the segmentation stage, we apply U-net model to get cropped tumor images. We input the original CT image to U-net and let it fit its corresponding mask. The architecture of U-Net is displayed in Fig. 2. We choose Dice loss (Eq. 2) here [16] to measure the contact ratio between predicted mask (P M ) and ground truth (GT ): where, p i is the pixel value in P M , t i is the pixel value in GT , N represents total pixels number of an image. And to avoid division by zero, a small constant is added to denominator.

Correlation between Gene and Image
To build the correlation between genes expression data and image features of tumor region in CT image series, we establish a conditional autoencoder to map the distinct source data for keeping the feature extraction ability without missing match. In our previous works, we find that there are two key problems: 1) The dimension of gene data can greatly influent the training effect so that leading to the model collapse.
2) The basic autoencoder can extract multiple level image features. Unfortunately, these features are only obtained without restriction by domain knowledge, so this usually leads to poor correlation between the gene and image feature. In this part, we introduce a gene encoder method for the dimension reduction and a conditional autoencoder for knowledge related feature extraction. Gene Encoder: The huge dimension of genomic data array is obviously an obstacle to array prediction, which is usually up to 10k. To handle this problem, like the operation in [40] for encoding word vectors, we introduced an encoder for gene data dimension reduction as follows.
Let gene data be a matrix g mn , m is the number of subjects and n is the length of gene array. We use the matrix w ns (s is objective dimension and is much smaller than n) to multiply g mn for encoding gene array. The equation is y = gw. Then we normalize the encoded gene array y to [0, 1] according different gene. Details of the gene encoder is listed in Fig. 4.
Conditional Autoencoder: When the output of encoder E is fitted to encoded gene y, the weights of E will be tended to transform image to y, but the image information will gradually vanish in this process. To avoid this situation, we add 1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  56  57  58  59  60  61  62  63  64  65 the decoder D after the E to fit original image X so that every layer of E can keep image information as much as possible. This process can be formulized as: whereỹ is the output of gene encoder, θ e is weight of E, θ d is weight of D and f eat is middle output of E. Thus, the conditional autoencoder can extract image features which have high correlation with gene data, and simultaneously keep images information. We apply this architecture to map gene and TR. While training the conditional autoencoder, we use perceptual loss [31] and L 1,2 loss to fitX to X: where φ is the pre-trained VGG16 CNN model trained on ImageNet, j is the j th layer in φ, here we use the layer named 'block5 conv3', and C j H j W j is feature map of j th layer. L 1,2 loss combines the L 1 loss and L 2 loss: We use this L 1,2 loss to fitỸ to Y : So the final loss function is L CA :

GCVAE-GAN for Pathological Visualization
Traditionally, radiogenomic methods can give a mathematic demonstration on the correlation between gene and images, but most of these methods merely provide the relationship numerically without any visual demonstration, and the visual results are greatly demanded for prognosis inferences.
The CVAE-GAN [35] encodes images and categories to Gaussian distribution, then decodes them to original images under the condition of classifier and discriminator. It can generate new images intra category by using noises due to there are many images belong to same category. But for radiogenomics, each tumor region CT image only corresponds to one gene array, which leads to the collapse of generation intra category. Focusing on the problem of data imbalance, we introduce a genomic conditional variational autoencoder GAN (GCAVE-GAN).
Specifically, we change the CVAE-GAN to fitting the tumor images as follows: 1) To solve data limitation, we let each subject belong to different class, and each class contains several different slices of tumor. So our model can generate different TRs by using the interpolation between genes instead of different noises. 2) We apply 1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  56  57  58  59  60  61  62  63  64  65 the encoded gene data as categories of corresponding TRs. For magnifying distance of encoded gene data in gene space, we set the threshold value to 0.5 and replace the values of gene with 0 and 1. Thus, gene data will perform as multi-class labels which control the generation process of model. See Fig. 4, our model is composed of four parts: 1) An encoder E which projects TRs to latent space z, 2) A generator G which transforms the latent vector z to TR, 3) A discriminator D which judges whether an image is real or fake TR, 4) A classifier C which projects TR to gene space.
To train this model well, there are six loss functions to support whole process [35]. To make the D have strong judgment ability to distinguish real TRs and synthesized TRs, D must minimize the loss function: To address the problem that the gradient of G is unstable while training, we use the mean feature matching objective proposed in [35]: which acts as perceptual loss, where f D (x) denote features on an first fully connection (FC) layer of the discriminator. It will accelerate the converging speed and improve stability of whole model.
Because the encoded gene data c show the tendency whether a group of genes is overexpressed or underexpressed, we don't fit it using logarithmic loss function. Instead, we apply the mean square loss on it.
To make the TR generated by G is belong to corresponding category c, G need to minimize the function The mean feature matching objective is also used here for the same purpose of L GD .
For the latent vector in each sample, we adapt the same strategy as [36]. The encoder outputs the mean µ and covariance of latent vector. Then use Kullback-Leibler (KL) loss to let distribution of latent space close to gaussian distribution. The KL loss can be formulated as Then latent vector z can be sampled as z = µ + r exp( ), where r N (0, I) is a random vector and represents the element-wise multiplication [35]. After z is inputted to generator G, we can get the synthetic TRx. We combine the L 2   1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  56  57  58  59  60  61  62  63  64  65 reconstruction loss, perceptual loss and pair-wise feature matching loss to minimize the difference between x andx, where L φ,j f eat is defined as Equation 5, f D and f C are the features of first FC layer of discriminator D and classifier C, respectively.
The final objective function is expressed as: where the formulations of each term are defined in Equation 11∼16. This model is trained as the original CVAE-GAN and the parameters are λ 1 = 3, λ 2 = 1, λ 3 = 10 −3 and λ 4 = 10 −3 .

Experiments and Results
In this section, we first conduct three experiments to separately demonstrate the results of models' output in segmentation, correlation and generation stages. Then, we prove the whole framework is comparable with traditional radiogenomics research on statistical analysis results. Finally, we show the generation results on the influence of gene change, and prove that the variance of gene expression will affect tumor states in radiology appearance.

Datasets
NSCLC Radiogenomics: Non-Small Cell Lung Cancer (NSCLC) is from Cancer Archive [30], where cancerous regions are relabeled as mask by doctors. In this paper, we choose NSCLC to construct our framework and evaluate the linkage among multi-source data. This dataset contains 211 subjects from a NSCLC cohort. Each subject in this dataset includes a CT series which contains hundreds of dicom images. And there are corresponding annotations of each tumor in each subject. Due to these annotations are just a coordinate of tumor, we invite doctors to relabel tumor as mask for us.

ROIs Extraction from CT Series
We separated the NSCLC dataset into training and testing set. The training set has 50 subjects up to 15000 CT images, which contains only image data for TRs segmentation [14]. The testing set has 161 subjects, each of them composed of CT image and its gene expression data. And the these 161 subjects was employed for the gene and image correlation latter. We choose U-net with Dice loss as the detector to recognize the position of tumor. We train it for just 500 epochs using Adam optimizer [32] with a learning rate 0.001. Then it is applied for testing set to detect position of tumor. TR is cropped at a size of 128×128. The cropped TRs and their corresponding gene data then form the dataset of correlation stage. Fig. 5 shows the predicted mask of UNet. The pixels belong to tumor are colored in blue which also show the tumor location. These masks cover main part of tumor   1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  56  57  58  59  60  61  62  63  64  65 except some irregular edges because of the less obvious edge of tumor, which is precise enough to the detection task. The size and location of tumor can be obtained at the same time. When we crop CT images, some surrounding information, such as bronchial tube and lung edge, which is crucial to mapping stage, are contained for providing more useful information.

Correlation between TRs and Genomic Data
In this part, we preprocess the gene data to remove the subjects whose gene have not available values. Finally, we get a dataset which has 113 subjects, and each of them contains 6 TRs of one tumor. We randomly sample 90 subjects as training set and 23 subjects as testing set, the ratio of testing set to training set is 1:4. The gene data are encoded using the gene encoder in Section . Then the conditional autoencoder can be trained as a multi-task network. Data augmentation is also employed in training process including horizontal/vertical flip and rotation. We train it using Adam optimizer [32] with a learning rate 0.001 until mean square error has decreased to 0.002 on testing set.
The images that conditional autoencoder outputs are shown in Fig. 6. The similarity between autoencoder's outputs and original TRs indicates model has already extracted images features sufficiently and won't be effected by geometric transformation. Meanwhile, it keeps the correlation between gene and image.
To prove the effectiveness of image features extracted from TRs, we apply hierarchical clustering on gene data and image features. We apply ResNet50 as the encoder of autoencoder. The features of different levels are extracted from last layer of each residual block. Then we apply Locally Linear Embedding (LLE) on them to reduce the dimension of features. Finally, we get four groups of features from different levels for each subject. The numbers of features are 64, 256, 512 and 1024, separately. Then we cluster image features among the samples to prove the correlation between features and disease (shown in Fig. 8). Like the work in [14], we also perform Chi-squared test to verify significant association of the DL features with prognosis data (shown in Table 1).
Also, to prove the effectiveness of condition, we remove it to do ablation experiment. We train the basic autoencoder without condition using the same process and extract features from same position of encoder. Then we calculate the correlation matrix of features (with and without condition) and gene data. Finally, we use box figure (Fig. 7) to show the difference between similarity with and without condition in each level.
To further prove the correlation between gene and image, we do gene-set enrichment analysis using GSEA software [37]. First, the Pearson correlation coefficient matrix c is calculated by features f and gene g. In each level, 10 features are selected according to s, the sum of all correlation coefficient's absolute value of each feature. Then genes will be sorted according to s to form pre-ranked files l, which are input of GSEA pre-ranked analysis. GSEA will find whether there is a gene set with known function appears on the top or bottom of l, which indicates features have positive or negative effect on that gene set. Finally, we will get the false discovery rate (FDR) and normalized enrichment score (NES) of each found gene set. The whole analysis process is expressed as Algorithm 1 and the correlation result is showed in Fig. 9.  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  56  57  58  59  60  61  62  63  64  65 In this part, we construct a linkage between DL features and genomic data from image to gene and prognosis domain, the results demonstrate this linkage has the same effectiveness as it in traditional radiology features, in addition, our proposed linkage can give an multilevel results unlike the handcrafted features. DL features work as traditional radiology features in this stage and the results of GSEA and prognosis validation are convincing.

TRs Generation under Gene Guidance
In this part, gene data are transformed to binary arrays which are called as controllers. Our task is to visualize TR using the corresponding controller. Each controller is the sum of all gene expression values weighted by different coefficients. We train this model for 4k epochs and get the cogent results.
The encoder of model first encodes gene and TRs to latent space, which is very close to the gaussian space. To prove the continuity of model when it generates TR, we use interpolation between two controllers to generate new TR whose state falls in between two corresponding TRs. Unlike the linear interpolation, we switch the state of controller one by one as our "controller interpolation", because value of controller only presents gene is overexpressed or underexpressed (0 or 1). The TRs generated by this method are showed in Fig. 10. We can find that generated images are changing gradually and continuously with controllers switched.
The generator then concatenates gene and a noise z to form an array and raise its dimension for generation process. We use the t-distributed stochastic neighbor embedding (TSNE) to visualize these arrays in latent space, part of which is shown in Fig. 11. We can observe that tumors with similar characteristics are projected closely in the latent space. Tumors in Group 1 are close to right edge of right lung and tumors in Group 2 are close to front edge of lung. Given this result, we will find the common ground of these similar tumors during generation process and how these controllers affect generation of tumors.
We select several groups of tumors, mainly about location and size. The controllers with activated state (whose value=1) in every sample of a group are extracted as the characteristic controllers. To research location controller, as a control, a template tumor is generated using a random noise z and a zero-like array which indicates all controllers are frozen. Then we activate the location controllers (including front, behind, left and right side) and use them to generate tumors (shown in Fig. 12 row  1 ). An enhancement can be seen in corresponding location, which means tumor is closer to thorax. To further prove the ability of location controllers, we freeze them while generating tumors with corresponding lung edge (shown in Fig. 13 (a,b,c,d)). We can easily find that when the corresponding controllers are frozen, outline of lung losses on some level in the corresponding position, which means tumor is far away from thorax.
We apply similar method when research size controller. They are activated one by one while generating tumors (shown in Fig. 12 row 2). The template tumor is bigger and bigger with size controllers activated gradually. To further prove the ability of size controllers, we freeze them while generating big tumors and activate them while generating small tumors. The results are shown in Fig. 13 (e,f), which indicates strong correlation between these controllers and tumor size.
In this part, we also link the genomic data and corresponding TRs. From the direction of gene raw data to TR images, we build a GCVAE-GAN model for this transformation. Considering the generation ability of model and the limitation of data, results of interpolation experiment shows that the model is continuous (see Fig. 10). In addition, we also change the gene expression values manually in a reasonable way. We validate the proposed framework by using the direction from genomic data to TRs. By this way, we use interpolation method to generate new TRs and change gene expression values to observe the variation of TR. The results prove that there is strong linkage from genome data to CT images. [14] proved that the genomic data and image computational features have strongly correlation with each other. We also demonstrate this result by using the proposed framework. And it will provide a new noninvasive way of investigating phenotypic information as traditional radiogenomics.

Discussion
After bi-direction validation(gene to image and image to gene), we find that some gene sets have high NES while correlated with DL features, among which, most gene sets are associated with cancer directly or indirectly by regulating cell activity (see Table 2). We cluster the samples into several distinct groups by using the deep features. Putting aside its result verification, we perform a statistical analysis on these cluster results and the prognostic data (T-stage and histology). The p value of Chi-squared test is low (which is ≤ 0.05) between them, see Table 1. This results demonstrate that our proposed framework can fulfill the correlation task all by using multi-levels deep features, especially on the correlation of tumor types and stages, that is the same as traditional methods [14].
The previous work [38] has demonstrated that the location of tumor has strong correlation with specific metagenes. In this context, we extent these anatomy information by some parameters, such as tumor location with thorax and tumor size. In this way, we find these genes as controllers that can guide the tumor growth. Different controller will enhance pixel values of corresponding part, which implicates there are some relations between gene and tumors. Due to each controller is consist of all gene values with different weights, an activated controller indicates high expression value of genes with high weights.
We change the location controllers and find that tumor will be away from thorax when froze the corresponding controllers (see Fig. 13(a, b, c, d)). In another way, tumor size is an important information in its prognosis and it can reflect the tumor states and chemoradiotherapy. We investigated the tumor size related genes sets and find that these size controllers have strong correlations with tumor growth. In our framework, we attempt to transform the tumor anatomy structure size by related genes, and the results indicate that the size controllers can enlarge and shrink the tumors when activated or frozen without disturbing other part of tumor region (see Fig. 13 (e, f) ).
Benefiting by this integrated dataset which has all data needed in radiogenomic researches, including original CT images, genomic and prognosis data, we can find and study the linkage between image and gene. This framework provides a method to find specific gene-set which has high correlation with specific DL features and gene which has high correlation with needed trait of TRs.
Limitations: Because of limitations relevant to the dataset scale and unbalance issues, our visualization results show that there are some abnormal instances because of lung edge and other organ like heart. The model can't ignore this obvious trait to pay attention to others. There would be more controllers which control other traits of tumor with a large dataset which is classified by tumor location. But in our work, these controllers are not completely independent because of the coupling phenomenon in tumor generation. Different groups of controllers influent each other and will change other traits of tumor.
As for the proposed framework has showed its capability in mapping the genomic data and tumor images, we will focus on dataset enlargement and combining different dataset from different sources and data categories.

Algorithm 1 Correlation Analysis
Input: The gene matrix, g; The feature matrix, f Output: The cluster heatmap between features and patients, h; The p value of Chi-squared test on predicted clusters and pronositic data, p; The GSEA analysis result of each level features and gene, r For i = 1 to 4 do 1: //Preprocess the features; 2: f i ← Transpose the last dimension of f i to second and flaten the last two dimensions; 3: f i ← Apply LLE on f i for reducing last dimension to 1; 4: //Correlation between patients and features; 5: h i ← Apply hierarchical clustering on f i ; 6: p i ← Use chi-square test between predicted clusters and prognosis; 7: //Correlation between g i and f i ; 8: //Now shape of f i is (n samples, n features), shape of g i is (n genes, n samples; 9: c i ← g i f i ; 10: s i ← sum(c i ) on gene dimension; 11: g i ← sort(g i ) according to s i ; 12: sg i ← Select top ten genes in g i ; 13: sg i ← sort(sg i ) according to c i ; 14: l i ← Combine gene identify with sg i ; 15: r i ←Input l i to GSEA; return h, p, r; 1  2  3  4  5  6  7  8  9  10  11  12  13  14  15  16  17  18  19  20  21  22  23  24  25  26  27  28  29  30  31  32  33  34  35  36  37  38  39  40  41  42  43  44  45  46  47  48  49  50  51  52  53  54  55  56  57  58  59  60  61  62  63  64