Synthetic whole-slide image tile generation with gene expression profile-infused deep generative models

Article Synthetic whole-slide image tile generation with gene expression profile-infused deep generative models Graphical abstract Highlights dThe RNA-GAN model is proposed as a solution for the RNA- to-image synthesis problem dA VAE reduces RNA-seq data’s dimensionality, preserving tissue differences dRNA-GAN achieves superior image quality compared with a traditional GAN dExperts rate RNA-GAN synthetic tiles higher than GAN synthetic tiles Authors Francisco Carrillo-Perez, Marija Pizurica, Michael G. Ozawa, ..., Luis Javier Herrera, Jeanne Shen, Olivier Gevaert Correspondence [email protected] In brief Carrillo-Perez et al. report RNA-GAN, a method that enables multi-modal integration of RNA expression and imaging data to generate synthetic tissue tiles. The model, trained on expression data, produces high-quality lung and brain tissue images preferred by experts and demonstrates imputation capabilities. Carrillo-Perez et al., 2023, Cell Reports Methods 3, 100534 August 28, 2023 ª2023 The Author(s). https://doi.org/10.1016/j.crmeth.2023.100534 ll Article Synthetic whole-slide image tile generation with gene expression profile-infused deep generative models Francisco Carrillo-Perez, 1,2 Marija Pizurica, 1,3 Michael G. Ozawa, 4 Hannes Vogel, 4 Robert B. West, 4 Christina S. Kong, 4 Luis Javier Herrera, 2 Jeanne Shen, 4 and Olivier Gevaert 1,5,6, * 1 Stanford Center for Biomedical Informatics Research (BMIR), Stanford University, School of Medicine, 1265 Welch Road, Stanford, CA 94305-547, USA 2 Computer Engineering, Automatics and Robotics Department, University of Granada, C. Periodista Daniel Saucedo Aranda, s/n, Granada, 18014 Granada, Spain 3 Internet Technology and Data Science Lab (IDLab), Ghent University, Technologiepark-Zwijnaarde 126, Gent, 9052 Gent, Belgium 4 Department of Pathology, Stanford University School of Medicine, 300 Pasteur Dr, Palo Alto, CA 94304, USA 5 Department of Biomedical Data Science, Stanford University, School of Medicine, Medical School Office Building (MSOB), 1265 Welch Road, Stanford, CA 94305-547, USA 6 Lead contact *Correspondence: [email protected] https://doi.org/10.1016/j.crmeth.2023.100534 SUMMARY In this work, we propose an approach to generate whole-slide image (WSI) tiles by using deep generative models infused with matched gene expression profiles. First, we train a variational autoencoder (VAE) that learns a latent, lower-dimensional representation of multi-tissue gene expression profiles. Then, we use this representation to infuse generative adversarial networks (GANs) that generate lung and brain cortex tissue tiles, resulting in a new model that we call RNA-GAN. Tiles generated by RNA-GAN were preferred by expert pathologists compared with tiles generated using traditional GANs, and in addition, RNA-GAN needs fewer training epochs to generate high-quality tiles. Finally, RNA-GAN was able to generalize to gene expression profiles outside of the training set, showing imputation capabilities. A web-based quiz is available for users to play a game distinguishing real and synthetic tiles: https://rna-gan.stanford.edu/, and the code for RNA-GAN is available here: https://github.com/gevaertlab/RNA-GAN. INTRODUCTION Biomedical data have become increasingly multi-modal, which has allowed us to better capture the complexity of biological processes. In the multi-modal setting, several technologies are used to obtain data from the same patient, providing a richer representation of their biological status and disease state. In current clinical practice, often demographic, clinical, molecular, and imaging data are collected on patients. Making these data modalities available helps advance the goals of precision medicine. 1,2 For example, DNA and RNA sequencing are now widely used for the characterization of patients with cancer. 3,4 Somatic mutation and gene expression profiles can be used to improve diagnosis, define disease subtypes, and determine the treatment regimen for patients with cancer. 5,6 Similarly, in pathology, tissue slides are the cornerstone for a variety of tasks. This includes primary diagnosis based on visual examinations by pathologists as well as treatment recommendations based on insights revealed by, e.g., immunohistochemistry stains. 7 Specifically for oncology, tissue slides are a valuable resource to observe morphological MOTIVATION The acquisition of multi-modal biological data for the same sample, such as RNA sequencing and whole-slide imaging (WSI), has increased in recent years, enabling studying human biology from multiple angles. However, despite these emerging multi-modal efforts, for the majority of studies, only one modality is typically available, mostly due to financial or logistical constraints. Given these difficulties, cross-modal data imputation and cross-modal synthetic data generation are appealing as solutions for the multi-modal data scarcity problem. Currently, most studies focus on generating a single modality (e.g., WSI), without leveraging the information provided by additional data modalities (e.g., gene expression profiles). Cell Reports Methods 3, 100534, August 28, 2023 ª2023 The Author(s). 1 This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/). ll OPEN ACCESS and texture changes, which reflect the tumor and its microenvironment. 7–9 Since the digitization of tissue slides to whole-slide image (WSI) data (specifically hematoxylin and eosin [H&E]- stained images), they have become a key data source for training deep learning models in a wide range of clinically relevant endpoints. 10 The association of molecular and morphologic patterns is also gaining interest in the research community, especially in the area of computational pathology. 11 In particular, the relationship between genomic features and WSI features has recently been demonstrated, with several studies showing that these two modalities are complementary. For example, morphological features from WSI data have been shown to associate with genomic mutations, gene expression profiles, and methylation patterns. 5,12,13 The effect of these variations are visually examined and can be spotted in the tissue both by pathologists and deep learning models. 5,14 However, some genomic variations are utterly rare, 15 limiting the available data. This fact slows downs the creation of machine learning models for their detection using a routinely obtained data source, such as H&E images. Furthermore, the variety of morphologic patterns that can be found is limited to the small pool of patients, losing a global perspective of the disease. Moreover, in terms of machine learning models, studies have shown that the integration of both modalities leads to an improvement in the performance of diagnostic and prognostic tasks in cancer when sufficient data are available. 5,16–19 However, both modalities are not always available due to financial or logistical constraints. For example, the Genome Express Omnibus (GEO) database 20 has numerous RNA sequencing (RNA-seq) datasets available, but few datasets have the corresponding WSIs. Similarly, most medical centers have large archives of tissue slides but not yet the means to generate matched gene expression data. New multi-modal datasets are being created to deal with these issues, 21 yet the problem still occurs for most clinical datasets. Thus, opportunities for training models that require multi-modal data are missed, slowing down progress in advancing precision medicine. 22,23 Data scarcity is a concerning problem in the machine learning community, especially in the context of recent successes for non-medical applications where huge amounts of data are available. 24,25 Specifically in biomedical problems, large and diverse cohorts are necessary to develop accurate clinical decision support systems that depend on machine learning algorithms. 26 To overcome the scarcity of heterogeneous annotated data in realworld biomedical settings, synthetic data are increasingly being considered. 27 Generative models, which impute synthetic data that are indistinguishable from real data, potentially offer a solution to deal with this issue. Within generative models, generative adversarial networks (GANs) and variational autoencoders (VAEs) have been widely used for multiple data generation tasks and have obtained exceptional performances in previous studies. 28,29 In both cases, the models learn a latent space to draw samples that cannot be discerned from real data. VAEs learn a latent space by solving the task of accurately reconstructing the original data using an encoder and a decoder, 30 and GANs are unsupervised generative models based on a game theoretic scenario where two different networks compete against each other. 31 The synthetic data generated by these models can expand the diversity of samples that a model is trained on, potentially increasing the predictive performance and also improving the model generalization capabilities. It is important to emphasize that this comes at almost zero cost once the model is trained, contrary to generating new data. In addition, synthetic data have the advantage that, when they can serve as a faithful representation of real patient data, they can be easily shared without any regulatory hurdles for protected health information. Several studies have focused on the generation of single-mo- dality synthetic data for both RNA gene expression and WSI data. For example, the generation of gene expression data has been mainly in the context of data imputation and has been researched by leveraging the latent space of VAEs. Qiu et al. showed that b-VAEs, a special case of VAEs, can impute RNA- seq data. 32 Similarly, Way et al. proposed a VAE trained on pan-cancer TCGA data that is able to encode tissue characteristics in the latent space and also leverages biological signals. 33 Recently, Vinas et al. presented an adversarial methodology for the generation of synthetic gene expression profiles that closely resemble real profiles and capture biological information. 34 The generation of high-quality WSI tiles has also been researched in recent years given the success of GANs in generating natural images. 35,36 For example, Quiros et al. showed that GANs are able to capture morphological characteristics of cancer tissues, placing similar tissue tiles closer in the latent space, while generating high-quality tiles. 37,38 However, current research focuses on generating or imputing single modalities without leveraging the information provided from other data types. For non-medical applications, crossmodal data generation has made enormous progress thanks to the availability of large paired data, e.g., paired text and image data. Unsupervised learning methods such as GANs, transformers, 39 and diffusion models 40 have been developed to leverage the relationship between these two modalities, enabling the generation of images based on their textual description 41–43 or generating textual descriptions of given images. 44 While cross-modal generation has proven successful for natural images in non-medical applications, 41,42 the relation between WSIs and gene expression needs yet to be explored for crossmodal synthetic data generation. For this use case, we were inspired by the observation that the relation between textual descriptions and their corresponding images is similar to the relation between WSIs and genomic data, as they are describing the same phenomenon from two different perspectives. Therefore, in this work, we develop an architecture for cross-modal synthetic data generation using WSIs and genomic information for heterogeneous healthy tissues as a use case. Specifically, we explore the generation of WSI tiles using gene expression profiles of healthy lung and brain cortex tissue (Figure 1). First, we train a VAE that reduces the dimensionality of the RNA-seq data. Then, using the latent representation of the gene expression as input, we present a GAN-based architecture (named RNA-GAN) that generates image tiles for healthy lung and brain cortex tissues. Based on evaluations of blinded pathologists, we show that the quality of generated WSIs can 2Cell Reports Methods 3, 100534, August 28, 2023 Article ll OPEN ACCESS be significantly improved when the GAN is infused with gene expression data. RESULTS Ab-VAE model can build a representative latent space that discriminates between healthy tissues As a first step, we aimed to create an accurate, distinguishable latent representation of healthy multi-tissue gene expression using a b-VAE architecture (Figure 1A). The goal was to reduce the dimensionality of the gene expression profile while maintaining the differences among the tissues. To do so, we use the traditional approach for training a b-VAE, i.e., reconstructing the input from the latent space (see STAR Methods). The b-VAE model was able to accurately reconstruct the gene expression by forwarding the latent representation through the decoder and obtaining a mean absolute error percentage of 39% (root-mean- square error [RMSE] of 0.631) on the test set for multiple tissues. To verify that the latent representation learnt by the b-VAE accurately maps to the different tissues, the uniform manifold approximation and projection (UMAP) algorithm 45 was used to visualize the real gene expression data as well as reconstructions of latent representations on the test set. For lung and brain samples, two separated clusters can be distinguished, showing how the model is characterizing the two tissues in the latent space (Figure 2A, ‘‘real’’ versus ‘‘reconstruction’’). Figure 1. Model architecture for gene expression, WSIs, and combined data using VAE and GANs (A) b-VAE architecture for the generation of synthetic gene expression data. The model uses as input the expression of 19,198 genes. Both the encoder and the decoder are formed by two linear layers of 6,000 and 4,096, respectively. The latent mand svectors have a feature size of 2,048. (B) GAN architecture for generating tiles by sampling from a random normal distribution. The architecture chosen was a deep convolutional GAN (DCGAN), 45 using as input a feature vector of size 2,048. The final size of the tiles generated is 256 3 256, the same as the size of the real tiles. (C) RNA-GAN architecture where the latent representation of the gene expression is used for generating tiles. The gene expression profile of the patient is used in the b-VAE architecture to obtain the latent representation. Then, a feature vector is sampled from a scaled random normal distribution (values ranging between [0.3,0.3]) and added to the latent representation. A DCGAN is trained to use this vector as input and generate a 256 3256 sample. The discriminator receives synthetic and real samples of that size. To further validate the learned latent space, we tested what happens when interpolating data in it. By interpolating in the latent space, we should be able to ‘‘transform’’ a randomly drawn sample to a gene expression profile that looks like it originated from one of the tissues (i.e., synthetic gene expression generation). To do so, we need to calculate the cluster centroid vector over the real data latent representations of the desired tissue and add this centroid vector to randomly drawn samples from the b-VAE latent distribution. This procedure allows us to generate synthetic gene expression data that look like real brain or lung gene expression data. When projecting these synthetic samples in the UMAP space, they indeed fall in the same clusters as the original data (Figure 2A, ‘‘generated’’ versus ‘‘real’’). We can also perform other operations in the latent space. For example, we should be able to ‘‘shift’’ the gene expression from one tissue into what it would look like if it originated from another tissue. In this case, weneedtoaddthedifference vectors between the cluster centroids of the respective tissues to the latent representation of a given sample gene expression. For example, we can shift a real brain gene expression profile to a lung gene expression profile, and vice versa. Visualizing these new samples in the UMAP space verifies that these operations can indeed be successfully performed (Figure 2B). Next, the representation capabilities of the b-VAE can also be extended to multiple tissues, showing a diverse representation with well-differentiated clusters and maintaining the generative capabilities across the multiple tissues (Figure 2C). The RMSE value distribution of the tissues is presented in Figure S3. Cell Reports Methods 3, 100534, August 28, 2023 3 Article ll OPEN ACCESS GANs generate quality synthetic WSI tiles preserving real data distribution differences Next, we developed a traditional GAN model to generate synthetic WSI tiles for brain cortex and lung tissue. The model was able to generate good-quality images, preserving the morphological structures and showing few artifacts (Figures 3A and S1A). In some tiles, checkerboard artifacts are noticeable, which is a known problem in GANs. 46 Despite the artifacts, the main cell types can be observed in the tiles, such as epithelial, connective, and muscle tissues. In addition, there is a clear distinction between the tiles generated for the brain cortex and for the lung, preserving the characteristics of the corresponding real tiles. Specifically, the brain cortex tissue is grouped in a set of layers that form a homogeneous and continuous layer (the outer plexiform layer, outer granular layer, outer pyramidal cell, inner granular layer, inner pyramidal layer, and polymorphous layer). 47 These characteristics can be observed in the synthetic brain tiles, i.e., they appear more homogeneous and contain less white spaces compared with the synthetic lung tissue tiles. The synthetic lung tissue tiles also present the characteristics of real tiles, showing the terminal bronchioles, respiratory bronchioles, alveolar ducts, and alveolar sacs in some cases. To test if the generated tiles have the same distribution as the real ones, the feature vectors outputted from one of the last convolutional layers of an Inception V3 network pretrained on ImageNet were obtained for the 600 generated tiles. Then, these feature vectors were projected and visualized using the UMAP al- Figure 2. UMAP visualization of b-VAE embedding of multi-tissue expression profiles (A) UMAP visualization of the real and reconstructed gene expression profiles of lung and brain cortex healthy tissues. Generated gene expression profiles, by sampling from the latent space and interpolating to the respective tissue, are also plotted. (B) Shifting real gene expression profiles between the two tissues. The latent representation of all the available samples is obtained, and the difference vectors between the cluster centroids are computed. (C) UMAP visualization of real gene expression profiles of multiple tissues and generated one from brain cortex tissue. gorithm, showing a similar distribution between the tissues for both real (Figure 3B) and synthetic samples (Figure 3C). Using latent gene expression profiles as input on GANs improves synthetic H&E tile quality and reduces training time Next, we used latent gene expression profiles as input instead of a random normal distribution for a GAN model generating WSI tiles. The gene expression was first forwarded through the pretrained b-VAE to reduce the dimensionality and to encode it in the latent space. Then, that representation plus noise sampled from a narrowed random normal distribution (values between [0.3,0.3]) was used as input to the generator, which outputs the synthetic tile (Figure 1C). This model generates synthetic tiles with fewer artifacts and better quality of the morphological structures (Figures 4A and S1B). To demonstrate that the gene expression latent representation provides actual information to generate the tiles and that the model does not mainly focus on the random signal (as with conventional GANs), we also created a GAN that samples only from a scaled random normal distribution (values between [0.3,0.3]). This model was not able to produce quality samples of any tissue (Figure 4C), even when sufficient training time was used. Hence, the RNA-seq data in the RNA-GAN is the main signal that guides its generation process and shows that the latent RNA-seq distribution is enough to generate synthetic tiles. The scaled normal distribution only adds sufficient noise to vary the generated tiles. We also obtained the feature vectors from one of the last convolutional layers of the Inception V3 architecture pretrained on ImageNet to observe if the distribution of the synthetic tiles was similar to that from real patients. The differences between the tissues were preserved, as well the tissue inner-cluster distribution (Figure 4B). To test the generalization capabilities of the trained models, we also used as input external brain cortex and lung tissue RNA-seq data (GEO: 120795). The model was able to 4Cell Reports Methods 3, 100534, August 28, 2023 Article ll OPEN ACCESS successfully generate tissue samples with characteristics similar to those obtained with the training data (Figure 4D). We then tested whether a model trained on real data can distinguish the synthetic generated tiles from this GEO cohort. This model reached an accuracy of 80.5%, an F1 score of 79.7%, and an area under the curve (AUC) of 0.805, showing that a model trained on real tiles can accurately classify the synthetic tiles. Finally, we observed that the RNA expression-infused GAN model needed fewer training epochs compared with the regular GAN model (Figure 5). To show the usefulness of these synthetic tiles in a clinical application, we decided to use them paired with self-supervised learning to show the improvements that can be obtained in a classification task between glioblastoma (GBM) and lung adenocarcinoma (LUAD) tiles. We firstly generated 10,000 tiles (5,000 per class) and trained a ResNet-50 architecture in a self-supervised learning (SSL) manner using the SimCLR framework. 48 SimCLR is a contrastive learning method that maximizes the agreement between two different augmented versions of the same image, learning a relevant feature representation of the image. Once the pretraining has been performed, the learned weights are used as the initialization weights for the downstream task and compared with training the model from scratch. Then, we took 1,000 real tiles from patients with GBM and LUAD from The Cancer Genome Atlas (TCGA) project, respectively, and performed a 5-foldcross-valida- tion in two settings: training from scratch, and using the SSL weights as initialization. The details of the SSL training are presented in the STAR Methods section. The model using SSL weights outperformed the one training from scratch in both terms of accuracy and F1 score. The model using SSL weights obtained a mean accuracy of 85.50% ±1.54% and an F1 score of 85.49% ± 1.54%, while the model trained from scratch obtained a mean accuracy of 76.71% ±4.51% and an F1 score of 76.67% ±4.52%. Both the accuracy and F1 score obtained across the splits and the confusion matrices are presented in Figure S4. Expert evaluation of synthetic tiles Next, we asked a panel of board-certified anatomic pathologists with different subspecialty expertise (M.G.O., H.V., R.B.W., C.S.K., Matt Van den Rijn, J.S.) to rate the quality of brain cortex and lung cortex tiles (N = 18) on a scale from 1 (worst) to 5 (best). These pathologists were not informed about the presence of synthetic data in the examples to prevent any bias when examining the tiles. Instead, they were told that the tiles presented were going to be used to train machine learning models and that we wanted to score their quality of that task. The set of questions asked and the procedure are further described in the STAR Methods section. Figure 3. A GAN generates realistic lung and brain cortex tiles maintaining the distribution of the real tiles (A) Tiles generated by the GAN model for brain tissue on the top and for lung tissue on the bottom. (B) UMAP representation of the real patients in the lung and brain cortex dataset. (C) UMAP representation of generated tiles using the GAN model. 600 tiles are generated per patient and then used to compute the feature vectors and the UMAP visualization. Cell Reports Methods 3, 100534, August 28, 2023 5 Article ll OPEN ACCESS The pathologists’ evaluations of the morphological structures resulted in a mean score of 3.55 ±0.95 for real brain, 2.88 ± 0.62 for GAN brain, and 2.94 ±0.64 for RNA-GAN brain. For the lung tissue, the mean score for the real samples was 2.26 ± 1.14, 1 ±0.55 for the GAN lung, and 1.73 ±0.79 for the RNAGAN lung. Hence, the pathologists rated the real samples as best quality, with second-best ratings for the samples from the RNA-GAN and the worst ratings for those from the conventional GAN. In the case of the RNA-GAN lung generated tiles, one of the pathologists scored the generated tiles higher than the real tiles in terms of quality (mean value of 1.2 compared with the score of real tiles, 0.8). Nevertheless, the rest of the pathologists scored the synthetic tiles lower than the real ones, but this difference was never greater than one point (real lung scores: [4,1.4,2,0.8,2.4,3]; RNA-GAN lung scores: [3,1,1.2,1.2,1.6,2.4]). The difference was greater in the GAN lung generated tiles, where the scores obtained were (1,0.6,1.2,0.6,0.6,2). This was not the case for the brain tissue, where pathologists scored similarly GAN and RNAGAN tiles, even though they scored slightly better the RNAGAN-generated tiles (real brain tile scores: [2.3,3.3,3,4.3,3.3,5]; RNA-GAN brain tiles scores: [2,3,3,3,2.7,4]; GAN brain tiles scores: [2.3,2.7,3,3,2.3,4]). The preference in the pathologists’ evaluations for the RNAGAN lung over GAN synthetic tiles is statistically significant (p = 0.025), while there was no statistically significant difference in ratings between the real and RNA-GAN lung tiles scores Figure 4. A gene expression-infused GAN improves tile quality (A) Tiles generated using the RNA-GAN model for lung and brain cortex healthy tissues. (B) UMAP visualization of the patients by generating tiles using their gene expression. The model preserves the distribution differences between the two tissues. (C) Generated tiles of model trained using only random Gaussian data on a small range ([0.3,0.3]) does not generate high-quality tiles, showing that the gene expression distribution is essential for synthetic tile generation. (D) Brain cortex and lung tissue tiles generated using an external dataset (GEO: GSE120795), showing the generalization capabilities of the model. 6Cell Reports Methods 3, 100534, August 28, 2023 Article ll OPEN ACCESS (p = 0.052). On the contrary, generated brain tissue tiles were significantly different from real tiles both for GAN (p = 0.043) and RNA-GAN (p = 0.038), and no significant difference in ratings was found between RNA-GAN and GAN (p = 0.305). However, a bigger mean evaluation score difference was found between real and GAN tiles than between real and RNA-GAN tiles, again confirming that the quality of RNA-GAN synthetic tiles is closer to the quality of real tiles (Figure 6A). In addition, the mean difference in the evaluation between GAN and RNA-GAN tiles was bigger than zero, showing the preference of pathologists for the RNA-GAN tiles over GAN tiles (Figure 6B). Finally, pathologists detected the tissue of origin of the tiles with a 100% ±0% accuracy for both real and RNA-GAN tiles, while this drops to 74.98% ±20.40% accuracy for the GAN tiles. DISCUSSION Previously, in biomedical problems, imputed samples and synthetic data were generated in an isolated way, not using the information provided by other modalities. In contrast, in this work, we studied the generation of WSI tiles of lung and brain cortex tissues by leveraging the corresponding gene expression profiles. This idea was inspired by recent advances in cross-modal data generation for non-medical images, including the generation or modification of images based on text prompts (e.g., DALL-E 2, 42 Imagen, 48 Parti, 49 and related models 41,42,44 ). Here, we extrapolated this idea to biomedical data by treating RNA-seq data as prior text to contextualize image generation. 50,51 This differs from other approaches presented in literature, where cost- and time-expensive transcriptomic profile levels were predicted Figure 5. A gene expression profile-infused GAN converges faster: Brain cortex and lung tissue tiles generated at the same epoch during training for the model with and without gene expression profiles The visualized epoch is the last epoch of training for the models using RNA-seq data. (A) Brain cortex generation at training epoch 24 for GAN and RNA-GAN models, with similar performance and quality between the generated tiles; however, less diversity is obtained when not using gene expression profiles. (B) Lung tissue generation at training epoch 11 for both the GAN and RNA-GAN models. A comparison of both models shows noticeable differences in the quality of the generated tiles. The model using gene expression profiles outputs better morphological features and less artifacts and has a higher overall quality. from an easy- and cheap-to-obtain paired H&E sample. 12 However, our approach can be situated in the digital biomarker study area and serve as a general purpose framework for other modalities. In addition, there is a growing archive of RNA-seq samples for which the H&E tissue samples have not been digitized missing that modality, and there are existing databases with publicly available gene expression data. Our proposed methodology can serve as an imputation method for generating the paired tissue slide from the RNA-seq profile of the patient. Several studies have been reported using GTEx data where gene expression variations have been studied across tissues and characterized several genetic variations. 52–55 These studies investigate how genetic variants affect gene regulation and complex traits in humans at the transcriptomic level, as well as the genetic bases of diseases. Thus, our proposed model can provide a solution to study those variations on tissues where only RNA-seq is available and model the effect that they might have on its morphology. Moreover, even though we are using healthy tiles as a use case to show our methodology capabilities, it can be further applied to other applications. The use of our method is potentially useful for rare diseases. In this context, synthetic data generation can be used as a tool to increase the availability of samples in scarce settings and can be used to study the variations between healthy and diseased H&E tiles. It is widely known that scarce data are one of the bottlenecks in the development of deep learning techniques 56 in medicine, and our proposed methodology can serve as a solution in this context. By pretraining models using recently developed SSL techniques, 57 global features can be learned with the synthetically generated data and then fine-tuned on the real samples. In addition, data quality has been presented in the literature as essential to obtain powerful machine learning models, 58 and in this work, we have shown how RNA-GAN- generated tiles have higher quality both in terms of pathologists’ evaluations and image quality metrics (e.g., Fre ´chet inception Cell Reports Methods 3, 100534, August 28, 2023 7 Article ll OPEN ACCESS distance [FID]) than those generated by a traditional GAN. In addition, the proposed generative framework could be further adapted to other imaging modalities including gene expression and computed tomography scans. The use of deep learning techniques is revolutionizing medicine, allowing the creation of systems that can accurately detect genetic variants or perform prognosis prediction. Nevertheless, for training these models, enough data are required, which can be especially challenging for rare mutation variants or rare cancer types (e.g., pediatric cancers). 56 Synthetic data generation can be used as a tool to increase the availability of samples in scarce settings, augmenting the number of samples. Using the generated data and SSL techniques, models can learn global features that can be further fine-tuned in downstream tasks with only a few samples. 57 However, with traditional generative approaches, it is not possible to regulate the effect of certain factors in tissue generation, apart from performing label guidance. This limits its potential in cases where we want to study the effect of an up- or downregulated gene in tissue morphology, expanding the diversity of samples. In this work, we proposed a methodology, RNA-GAN, that can be used under this assumption when trained on the required data, creating in silico simulations of specific genomic variations and studying their effect on the tissue. An example would be the generation of prostate cancer tiles where the TP53 mutation can be found, helping to increase the performance of detection models. 14 Moreover, the effect of specific mutations in a specific gene on morphologic characteristics of tissue has been studied but not from the expression of multiple genes. 59,60 Given that the latent representation of the gene expression is being used, the effect of thousands of genes can be reflected in the generated tissue tile, not limited to specific expression values. From a machine learning perspective, since biomedical data are becoming increasingly multi-modal, there is a growing interest in developing multi-modal predictive models for advancing the goals of precision medicine. However, obtaining multi-modal biological data is a slow and costly process. Multimodal data are sometimes available for widely studied diseases, but this is not common and definitely not the case for pediatric and rare diseases. In addition, not all medical facilities have the required expertise or instruments to collect each data modality for a patient. Therefore, cross-modal data imputation 32,61,62 and generation of cross-modal synthetic data is a promising approach to complete datasets by leveraging the potential of deep learning models. 27 As a first step in our work, we developed a b-VAE model that reduces the dimensionality of gene expression data. We have shown that the b-VAE model is able to capture the latent representation of multiple human tissues and that it obtains a representative latent feature vector that accurately distinguishes between the tissues (Figure 2). To test the accuracy, practicality, and capabilities of the learned latent space, we performed certain sanity checks. For example, we tested the ability to generate synthetic gene expression profiles and to interpolate samples between classes. Our analysis shows that the b-VAE model indeed allows us to obtain a compact, accurate representation of tissue gene expression profiles. We later use the encoder part of the b-VAE to reduce the dimensionality of gene expression profiles, which in turn will guide the generation of the WSI tiles. Note that this compact representation could also be used for downstream tasks including prognosis or treatment outcome prediction, but this analysis was out of the scope of this work. Next, we trained a traditional GAN model on brain cortex and lung tissue data to generate WSI tiles. We explored the possibility of using the same architecture to generate both tissues, but the model collapsed and only generated brain cortex tiles. We hypothesize that this could be due to the homogeneity of brain tissue compared with lung tissue, making it easier for the model to generate brain cortex tiles and reduce the loss. This can also be observed in the number of necessary epochs to generate quality tiles for each tissue. While the model only needs 34 epochs for training convergence for brain cortex tissue, it requires 91 epochs for lung tissue. Clearly, it is more difficult to generate the lung tiles, probably because lung tissues are more heterogeneous compared with the brain cortex. Figure 6. Expert evaluation of synthetic slides (A) Difference in morphological structure quality of synthetic (generated by GAN and RNA-GAN) and real tissues based on the pathologists’ evaluations. The difference between real tiles and generated tiles was bigger for GAN than for RNA-GAN. (B) Difference in morphological structure quality between the synthetic generated tiles by the GAN and RNA-GAN based on the pathologists’ evaluations. Pathologists evaluated the tiles generated using RNA-GAN better compared with only GAN. 8Cell Reports Methods 3, 100534, August 28, 2023 Article ll OPEN ACCESS of a certain dimension are taken from the tissue, and they are used to train the models, which is consistent with related work in state- of-the-art WSI processing. 5,71,72 In our work we took nonoverlapping tiles of 256 3256 pixels. Firstly, a mask of the tissue in the higher resolution of the SVS file was obtained using the Otsu threshold method. 73 Tiles containing more then 60% of the background and with low-contrast were discarded. A maximum of 4,000 tiles were taken from each slide. For the preprocessing of the images we relied on the python package openslide, 74 that allows us to efficiently work with WSI images. The tiles were saved in an LMDB database using as index the number of the tile. This approach enables to reduce the number of generated files, and structure the tiles in an organized way for a faster reading while training. bVAE architecture for synthetic gene expression generation and experiments We chose the bVAE model for the generation of synthetic gene expression data. 75 The bVAE model is an extension of the VAE where a bparameter is introduced in the loss function. The original auto-encoder is formed by two networks, the encoder and the decoder. The encoder encodes the input into a lower dimensionality representation, and then it is used to reconstruct the input using the decoder, by learning the function h q (x)zxbeing qthe parameters of the neural network. To learn this function, we want to minimize the reconstruction error between the input and the output. The most common loss function is the root mean squared error (RMSE). However, for the VAE we want to learn a probability distribution of the latent space, which allows us to later sample from it to generate new samples. The assumption of the VAE is that the distribution of the data x,P(x) is related to the distribution of the latent variable z, P(z). The loss function of the VAE, which is the negative log likelihood with a regularizer is as follows: Liðq;4Þ=EzqqðzjxiÞ½log p4ðxijzÞ +KLðqqðzjxiÞkpðzÞÞ (Equation 1) where the first term is the reconstruction loss and the second term is the Kullback-Leibler (KL) divergence between the encoder’s distribution q q (z|x) and p(z) which is defined as the standard normal distribution p(z)=N(0,1). For the bVAE we introduce the parameter b, which controls the effect of the KL divergence for the total loss: Liðq;4Þ=EzqqðzjxiÞ½log p4ðxijzÞ +b3KLðqqðzjxiÞkpðzÞÞ (Equation 2) If b= 1, we have the standard loss of the VAE. If b= 0, we would only focus on the reconstruction loss, approximating the model to a normal autoencoder. For the rest of the values, we are regularizing the effect of the KL divergence on the training of the model, making the latent space smoother and more disentangled. 75 For the final architecture, we empirically determined to use two hidden layers of 6,000 and 4,096 neurons each for both the encoder and the decoder, and a size of 2,048 for the latent dimension. Given that we were going to use the latent representation for the generation of the tiles, we followed the same dimensionality as the output of the convolutional layers of state-of-the-art convolutional neural networks. 76 We used batch norm between the layers and the LeakyReLU as the activation function. A b= 0.005 was used in the loss function. We used the Adam optimizer for the training with learning rate equal to 5 310 5 , along with a warm-up and a cosine learning rate scheduler. We trained the model for 250 epochs with early stopping based on the validation set loss, and a batch size of 128. A schema of the architecture is presented in Figure 1A. We divided the dataset in 60-20-20% training, validation and test stratified splits. We trained two different models, one for brain cortex and lung tissue data, and the other with all the tissues described in previous subsections (lung, brain cortex, stomach, pancreas, and liver). GAN architecture for synthetic WSI tiles generation and experiments GANs have been succesfully used for generating high-fidelity images. In this work we use the Deep Convolutional GAN architecture presented by Radford et al. 77 On early experiments we used the minmax loss function described on the original Godfellow et al. 31 work. However, this loss function led to a lack of diversity in the generation of the samples, and a diminished quality. Therefore, we decided to use the Wasserstein loss introduced by Arjovsky et al., 78 also adding the gradient penalty proposed by Gulrajani et al. 79 In this case the discriminator (or critic, as called in the paper) does not classify between real and synthetic samples, but for each sample it outputs a number. The discriminator training just tries to make the output bigger for real samples and smaller for synthetic samples. This simplifies the loss function of both networks, where the discriminator tries to maximize the difference between its output on real instances and its output on synthetic instances as follows: LD=E~ xPg½Dð~ xÞ  ExPr½DðxÞ +lE^ xP^ xhkV^ xDðb xÞk212i(Equation 3) and the generator tries to maximize the discriminator’s output for its synthetic instances as follows: LG=E~ xPg½DðGð~ xÞÞ (Equation 4) We trained two DCGANs, one per tissue, by sampling from a normal random distribution (scheme depicted in Figure 1B). We samples different number of tiles per image for the training of the network, finally selecting 600 tiles per image because the quality of the image and the artifacts were highly improved by augmenting the number of tiles. We used the Adam optimizer for both the generator and the discriminator, with a learning rate equal to 1 310 3 for the generator, a learning rate equal to 4 310 3 for the discriminator and betas values (0.5,0.999) in both cases. Data augmentations such as color jitter and random vertical and horizontal flips were used Cell Reports Methods 3, 100534, August 28, 2023 e2 Article ll OPEN ACCESS during training. The brain tissue GAN was trained during 39 epochs while the lung tissue GAN was trained during 91 epochs. For the training of the GANs, the Python package Torchgan was used. 80 RNA-GAN architecture for synthetic WSI tiles generation and experiments After the successful generation of the tiles using a traditional GAN approach, we explored the generation of synthetic tiles by using the gene expression profile of the patient. We combined the pretrained bVAE with the DCGAN architecture, using the encoding in the latent space as the input for training the generator. To generate different tiles from the same gene expression profile, we sample a noise vector from a narrowed random normal distribution (values ranging between [0.3,0.3]) and add it to the latent encoding. Therefore, the input to the generator would be: ~ x=qqðzjxÞ+Nð0;1Þ(Equation 5) We trained two DCGANs, one per tissue, and the pipeline is depicted in Figure 1 C). We finally selected 600 tiles per image to train the generator. We used the Adam optimizer for both the generator and the discriminator, with a learning rate equal to 1 310 3 for the generator, a learning rate equal to 4 310 3 for the discriminator and betas values (0.5,0.999) in both cases. Data augmentations such as color jitter and random vertical and horizontal flips were used during training. The brain tissue GAN was trained during 24 epochs while the lung tissue GAN was trained during 11 epochs. For the training of the GANs, the Python package Torchgan was used. 80 To validate the generalization capabilities of the trained model, the GEO series 120795 was used. It contains gene expression profiles from healthy tissues, where we took the expression of lung and brain cortex tissues. For obtaining machine learning performance metrics, one hundred images were generated per tissue. Then, a Resnet-18 was trained from scratch using 10 epochs and early stopping based on a 20% of data as validation set. A learning rate value of 3e 5 and AdamW optimizer were used. Finally, the model was tested on the GEO synthetically generated data, and accuracy, F1-Score and AUC was computed. To test the use of the synthetic tiles for pretraining machine learning models, we used SimCLR to train a Resnet-50 backbone on 10000 synthetically generated tiles (5000 per class). We trained the model for 100 epochs, using the SGD optimizer with a learning rate equal to 6 310 2 , a momemtum value of 0.9 and a weight decay equal to 5 310 4 , and we used a cosine annealing learning rate scheduler. For both training from scratch and fine-tuning the SSL weights, we used AdamW as the optimizer, with a weight decay of 0.01 and a learning rate value of 3 310 5 . Both models were trained for 40 epochs using a batch size of 4, and a stratified 5--Fold CV was used. We release an online quiz where users can try to distinguish between real and synthetic samples, obtaining a score on how well they performed. The quiz is available in the following URL: https://rna-gan.stanford.edu/. The code and checkpoints for the proposed models is available in this Github repository: https://github.com/gevaertlab/RNA-GAN. Expert pathologist evaluation form To evaluate the quality of the synthetic tiles, we presented a form to expert pathologists (Figure S2). The pathologists were not informed that some presented tiles were synthetic, to omit any kind of biases in the evaluation. Instead, we informed the pathologists that these tiles were going to be used to create machine learning classifiers, and we wanted to evaluate their quality for this task. Three questions were asked to the experts. 1. Is the tile from brain cortex or lung tissue? 2. Quality of the morphological structures: Being 1 very bad and 5 very good, how would you rate the morphological features present in the tile for an assessment of the tissue? 3. Do you find artifacts in the image? (e.g. image aberrations) (Yes/No) Ten tiles were randomly selected from 600 synthetically generated tiles per architecture. Then, tiles were visually inspected and the three with the highest quality were selected, leaving with a total of 18 tiles (between real, GAN generated, and RNA-GAN generated). This number of tiles was selected to enable clinicians to complete the form in less than 10 min. We used the same selection criteria for both cases (traditional GAN and RNA-GAN), and we manually selected them so both synthetic tiles could be reviewed in the best condition possible. Since synthetic tiles generated by the networks could present checkerboard artifacts (in less extend in RNA-GAN generated tiles), it was something that we wanted to omit therefore we used the following criteria. dNo checkerboard artifacts were present in the tile. dTissue characteristics were clearly visible in the tile. QUANTIFICATION AND STATISTICAL ANALYSIS Model performance was primarily assessed using expert pathologist’s evaluations, FID scores, Accuracy, and F1-Score. Further details can be found in the methods details and Figure Legends. e3 Cell Reports Methods 3, 100534, August 28, 2023 Article ll OPEN ACCESS