Interactive analysis and exploration of cell-cell communication in single-cell transcriptomics with InterCellar

Deciphering cell-cell communication is a key step in understanding the physiology and pathology of multicellular systems. Recent advances in single-cell transcriptomics have contributed to unraveling the cellular composition of tissues and enabled the development of computational algorithms to predict cellular communication mediated by ligand-receptor interactions. Despite the existence of various tools capable of inferring cell-cell interactions from single-cell RNA sequencing data, the analysis and interpretation of the biological signals often require deep computational expertise. Here we present InterCellar, an interactive platform empowering lab-scientists to analyze and explore predicted cell-cell communication without requiring programming skills. InterCellar guides the biological interpretation through customized analysis steps, several visualization options and the possibility to link biological pathways to ligand-receptor interactions. By analyzing COVID-19 and melanoma cell-cell interactions, we show that InterCellar resolves data-driven patterns of communication and highlights molecular signals through the integration of biological functions and pathways.

SingleCellSignalR 14 ), or CCI-data result of custom algorithms (e.g., Kumar et al. 15 ) that predict cell-cell interactions. The second step allows the user to focus the analysis on three biological domains, named universes: cell clusters, genes and biological functions. In each universe (cluster-verse, gene-verse and function-verse), different ltering and visualization options are provided, giving the opportunity to interactively subset the data and get qualitative and quantitative insights. Particularly, in the so-called function-verse, InterCellar provides multiple resources to perform functional annotation of enriched interactions. Lastly, InterCellar enables a focused inspection of cell-cell communication through the analysis of the so-called int-pair modules, which are de ned as groups of interactions with similar functional patterns. The user may select one cell cluster as viewpoint and one ow of communication among: (i) directed-outgoing, i.e., the viewpoint cluster sends ligands to other clusters; (ii) directedincoming, i.e., the viewpoint cluster expresses receptors; and (iii) undirected, as in the case of receptorreceptor pairs (see Methods). A subset of interactions is then used to de ne int-pair modules that can be visualized and analyzed together with their signi cant functions. This last step can be iterated for all cell clusters and ows of interest. Importantly, InterCellar ensures reproducibility of the analysis and facilitates collection and nalization of results by providing multiple download options for tables and gures. In the following, we present InterCellar's main features and demonstrate its general applicability on two datasets, concerning coronavirus disease 2019 (COVID-19) and melanoma.

Results
InterCellar highlights data-driven patterns of cell-cell communication in scRNA-seq data In order to demonstrate InterCellar's functionalities implemented in the three biological domains (socalled, universes), we considered a publicly available COVID-19 dataset from Chua et al. 6 , containing nasopharyngeal and bronchial samples from 19 patients and from ve healthy controls. In particular, we adopted the clinical classi cation of patients provided by the study, consisting of 8 moderate cases and 11 critical cases. Moreover, we retained the original cell type labeling performed by the authors and ran CellPhoneDB (v2.0) 11 to obtain predicted cell-cell interactions as input to InterCellar. Thus, three separate input datasets were generated, corresponding to control, moderate and critical cases. Cell types can be grouped in (a) epithelial cells, namely basal, ciliated, ciliated-differentiating (ciliated-diff), FOXN4+, ionocyte, IFNG responsive cell (IRC), secretory, secretory-differentiating (secretory-diff) and squamous and (b) immune cells, namely B cell, cytotoxic T cell (CTL), regulatory T cell (Treg), natural killer (NK), natural killer T cell (NKT), natural killer T cell-proliferating (NKT-p), neutrophil (Neu), plasmacytoid dendritic cell (pDC), monocyte-derived dendritic cell (moDC), mast cell (MC), resident macrophage (rMa), non-resident macrophage (nrMa) and monocyte-derived macrophage (MoMa).
As the rst step, we focused on the analysis of the cluster-verse and considered the total number of interactions per cell type while comparing critical to moderate cases (Fig. 2a). We observed rather small differences in numbers of interactions for cell clusters belonging to epithelial cells (with the exception of ionocytes), while immune cell types showed a higher variability. In particular, we could con rm the ndings of Lin et al. 7 , who described a higher number of interactions for macrophages (MoMa and nrMa) as well as a lower number of interactions for T cells (CTL and Treg), when comparing critical to moderate cases. Moreover, in the same comparison, comprehensive consideration of all cell types highlights a striking gain of interactions for mast cells (MC) as well as proliferating natural killer T cells (NKT-p). This novel nding holds true when examining the number of interactions between a certain cell type of interest and all other cell types ( Fig. 2b and Supplementary Fig. 1a,b). Speci cally, both immune cells (e.g., MoMa, nrMa, CTL and Treg) and epithelial cells (e.g., ciliated, ciliated-diff, secretory and secretory-diff) show a consistent pattern of communication in critical cases, characterized by a higher number of interactions occurring between all cell types and MC or NKT-p cells. Interestingly, MCs have recently been hypothesized to have a major role in driving hyperin ammation in severe cases of COVID-19, due to their dysfunctional phenotype related to mast cell activation syndrome 16,17 .

InterCellar unravels the composition of cellular communication by examining ligand-receptor pairs
Another biological domain relevant for cell-cell communication is represented by the proteins (and thus genes, in the case of scRNA-seq) that compose ligand-receptor pairs. InterCellar's gene-verse offers the possibility to investigate precisely which cluster-pairs communicate through which genes. We call "intpair/cluster-pair couplet" the occurrence of an int-pair in a certain cluster-pair. These couplets are visually represented by InterCellar in dot plots, where the analyst interactively chooses int-pairs and cluster-pairs to examine. Thus, we proceeded in the analysis of the COVID-19 datasets, by plotting the interaction score (calculated by CellPhoneDB, see Methods) of int-pairs occurring in the communication between secretory cells, secretory-diff and all other cell types. In particular, by looking at int-pairs composed of selected chemokine ligands mentioned in Chua et al. 6 (CXCL1, CXCL3, CXCL6, CXCL16, CXCL17), we could observe a clear enrichment of ligand-receptor pairs in the communication between secretory as well as secretory-diff cells and neutrophils (Fig. 3a). Speci cally, these int-pairs were detected at varying levels of expression in the two disease conditions, while they were completely absent in control samples. This nding con rms the hypothesis of the authors concerning neutrophil recruitment induced by secretory cells 6 . However, no relevant difference could be recognized between moderate and critical cases.
To further investigate the latter phenomenon, we made use of InterCellar's gene-verse functionalities by highlighting only int-pair/cluster-pair couplets that are unique to a certain phenotype (i.e., control, moderate or critical). We considered two groups of chemokine ligands, namely the CCand CXCsubfamilies and selected all possible int-pairs found in each phenotype. CXCL-pairs showed a predominant enrichment of interactions unique to critical cases, in both epithelial (47%) and immune cells (57%) ( Fig. 3b and Supplementary Fig. 2a Fig. 2a). For CCL-pairs, we focused on selected immune cell types, separated into three groups: macrophages (MoMa, nrMa and rMa), NK cells (NK, NKT and NKT-p) and T cells (Treg and CTL) (respectively Fig. 3c and Supplementary Fig. 2b,c). While for macrophages the proportion of unique interactions favors the critical phenotype (56%), an inverse tendency could be noticed for NK (41% critical) and T cells (20% critical). In particular, among macrophages, MoMa and nrMa displayed many chemokine interactions that were unique to critical cases. These interactions involved recipient cell types such as MC, moDC, neutrophils and NKT-p. Speci cally, MoMa communicates with these other cell types via CCL4L2 while nrMa via CCL7. At the same time, interactions of rMa via CCL3 in controls are absent in COVID-19 cases (Fig. 3c). On the contrary, NK and NKT as well as Treg and CTL showed an enrichment of interactions unique to moderate cases, directed towards secretory cells (among others) ( Supplementary   Fig. 2b,c). Altogether, these results are in line with two ndings by the authors: on the one hand, critical cases are characterized by a highly in ammatory pro le for MoMa and nrMa; on the other hand, a wellbalanced immune response distinguishes moderate cases, in which the communication between immune cells and epithelial cells underlie an effective response to the viral infection 6 .
InterCellar links biological functions and signaling pathways to cellular communication InterCellar's function-verse investigates the third biological domain of interest and offers the opportunity to annotate cell-cell communication with biological functions and pathways. We performed functional annotation for each COVID-19 dataset, using all functional databases provided by InterCellar (see Methods). Then, we examined functional terms of interest through a sunburst visualization, which combines all int-pair/cluster-pair couplets that have been annotated to that particular function. To investigate signals of a "cytokine storm" 18 in critical COVID-19 cases, as previously suggested by the analysis of chemokines, we selected the functional term "viral protection interaction with cytokine and cytokine receptor" (annotated from KEGG 19 ) and compared the sunburst plots of moderate and critical cases (Fig. 4). As a rst observation, we noticed that moderate cases had a higher total number of intpairs annotated to this function (681 pairs, compared to 553 in critical cases). However, when considering only unique int-pairs, critical cases had the highest number (51 vs. 46). While CTL and Treg are predominant in moderate cases (with the highest total int-pairs, respectively 104 and 101), these two cell types are replaced by MoMa and nrMa in critical cases, occupying third and fourth place in the ranking with almost half of their int-pairs (55 and 54, respectively). Interestingly, control cases are comparable to moderate ones with respect to the ranking of CTL, B cells and Tregs ( Supplementary Fig. 3). Furthermore, when examining the speci c int-pairs, CTL communicated in moderate cases with nrMa (as well as other cell types) through a cytokine-encoding ligand called LTA (lymphotoxin-alpha), which is completely absent in critical cases. On the other hand, MoMa and nrMa pro les seem to be fairly conserved in the two disease phenotypes, with only a few int-pairs unique to critical cases (e.g., IL18 & IL18 receptor for MoMa::CTL; CCL7 & CCR1, CXCL5 & CXCR1, CXCL5 & CXCR2 for nrMa::Neu). Overall, these results suggest once again a dysfunctional immune response in COVID-19 critical cases, where a balanced lymphocyte-mediated regulation seems to be overthrown by excessive activation of macrophages.
InterCellar guides a focused analysis of cellular communication through the de nition of functionally distinct modules The work ow of InterCellar concludes with the analysis of int-pair modules, de ned as groups of functionally similar interactions. In order to demonstrate this last step, as well as the general usability of InterCellar, we considered a second scRNA-seq dataset, composed of metastatic melanoma samples from Tirosh et al. 20 . This dataset includes both malignant cells and cells from the tumor microenvironment (TME) of 19 patients. We retained the cell type labels assigned by the authors, namely melanoma-malignant cells, T cells, B cells, macrophages (Macro), endothelial cells (Endo), cancerassociated broblasts (CAF) and natural killer (NK) cells. Once again, we ran CellPhoneDB (v2.0) 11 to obtain ligand-receptor interactions as input to InterCellar.
Since cancer cells have been described as actively recruiting and reprogramming normal cells of the tumor ecosystem 21 , we used InterCellar to perform an in-depth analysis of the interactions that characterize malignant cells in their communication with the TME. Thus, we chose the malignant cell cluster as viewpoint in the analysis and focused on directed-outgoing interactions as those int-pairs where the ligand is expressed (and sent) by malignant cells to receptors expressed on all other cell types.
Brie y, InterCellar's de nition of int-pair modules is based on the functional annotation performed in the function-verse: each module is composed of a subset of int-pairs that share common functional patterns.
These modules can be visualized in a two-dimensional plot as a UMAP 22 , which represents int-pairs clustered by functional similarity (Fig. 5a, see Methods). Moreover, for each module, InterCellar displays a circle plot to identify precisely which cluster-pairs and genes are involved (Fig. 5b-e). Hence, a total number of nine int-pair modules was de ned for the selected interactions and, for each module, we could associate functional terms found to be statistically signi cant (using a one-sided permutation test, see Methods). Interestingly, by examining the circle plots, the relevance of many int-pairs could be validated by further literature research. For example, we investigated the underlying interactions for modules #1, #3, #4 and #7. Int-pair module #1, characterized by the functional terms "extracellular matrix (ECM) organization" and "integrin signaling pathways", is composed of collagen-encoding genes expressed by malignant cells, which interact with integrin-complexes expressed by endothelial cells, CAF, B cells and T cells, as well as with malignant cells themselves. Notably, the importance of collagen-integrin interactions in promoting cell invasiveness has been already described in the literature by Zhou et al. 23 . EPH-Ephrin signaling, which distinguishes, among others, int-pair module #3, has been associated with tumor neovascularization 24 and is mainly promoted by interactions between malignant cells and endothelial cells or CAFs. Int-pair module #4 shows greater heterogeneity and comprises interactions associated with vascular endothelial growth factors (VEGFA and VEGFB), platelet-derived growth factors (PDGFA) as well as broblast growth factor (FGF5) that are enriched to a higher extent between malignant cells and CAFs, compared to the other cell types. Notably, the role of malignant cells in promoting CAFs proliferation has been previously investigated by Izar et al. 25 . However, the authors focused on the interaction between transforming growth factor-alpha (TGF-alpha) and epidermal growth factor receptor (EGFR), which was not found in this data. Finally, mechanisms of metastatic invasiveness in melanoma have been associated with the overexpression of annexin A1 in malignant cells 26 . Interactions involving this gene can be seen in the circle plot of int-pair module #7: ANXA1, sent by malignant cells, interacts with formyl peptide receptor (FPR) exclusively expressed by macrophages, suggesting a deleterious cross-talk between these two cell types.

Discussion
In the context of single-cell transcriptomics, the investigation of cellular cross-talk has been regarded as a valuable method for understanding biological mechanisms in health and disease 10 . Despite the existence of many algorithms for predicting the occurrence of cell-cell communication based on ligand-receptor interactions, further analysis and interpretation of the associated biological signals are not trivial due to the complexity of results and the lack of a standardized approach to perform such analyses 13 .
We developed InterCellar as an interactive analysis platform designed to guide scientists in the interpretation of complex cell-cell communication results obtained from scRNA-seq data.
By demonstrating InterCellar's functionalities on two different scRNA-seq datasets, we were able to obtain comprehensive insights on the modalities of cell-cell communication and the underlying ligand-receptor interactions, which could be further corroborated by literature review. Moreover, we have shown that InterCellar is a powerful tool for identifying and highlighting previously unknown molecular interactions inferred from single-cell transcriptomic data. For example, an in-depth analysis of cellular interactions in COVID-19 patients con rmed the primary role of immune response in discriminating between different courses of the disease 6 . In addition, the analysis suggested a novel, consistent pattern of communication for patients with critical COVID-19, involving cell types such as mast cells. These results may provide new exciting hypotheses to deepen the understanding of severe SARS-CoV-2 cases.
Implemented in R/Shiny, InterCellar is distributed as an open-source Bioconductor package, ensuring software robustness, maintenance and compatibility with multiple operating systems. With an existing installation of R and Bioconductor on the user's personal computer, installing the InterCellar package and launching the Shiny app can be accomplished by running two lines of code, which, together with further information, can be found in the online user guide (https://bioconductor.org/packages/InterCellar/). Notably, the local installation of InterCellar avoids issues related to data privacy or sharing of unpublished patient data, which is a typical problem for web-based platforms running on external servers. Moreover, each step of InterCellar's work ow is generally completed in a very short time frame (from seconds to few minutes), ensuring a truly interactive and performant analysis.
InterCellar directly builds upon the results of other existing tools developed to infer cell-cell communication, which can be freely chosen by the user. Alongside the automatic import from supported tools, InterCellar accepts custom data with speci c format and necessary information (see Methods).
Since prediction methods often rely on different reference databases, and build their results on diverse statistical and mathematical assumptions 10 , evaluating advantages and disadvantages of each method is critical for choosing the method that best ts to the data of the user. Lastly, we want to stress the importance of further experimental validation of the predicted interactions, which could be achieved, for example, by single-molecule uorescence in situ hybridization or measurement of co-occurrence through ow cytometry 10 .
To benchmark InterCellar's features, we identi ed 18 published tools from two recent review papers 10,13 as well as from manual search (as of April 2021, Supplementary Table 1). We collected methods allowing the analysis of cell-cell communication and excluded 13 methods that required programming skills (e.g., extensive coding in R/Python to perform the analysis) from a systematic comparison, thus focusing on the remaining web-based or standalone systems (Fig. 6). We found that one tool, CCCExplorer 27 , had comparable functionalities to InterCellar, providing a local, interactive analysis that focuses on enriched biological pathways. However, this software was developed for bulk RNA sequencing, therefore lacking a straightforward application in scRNA-seq. Three platforms, namely pyMINEr 28 , CellChat Explorer 12 and Cellinker 29 lack a proper interactive analysis, the possibility of customizing the data and download options; therefore, they could be regarded mostly as explorative tools providing static results. Moreover, CellChat Explorer and Cellinker run remotely on a server, thus requiring data sharing and possibly leading to privacy issues. The authors of CellChat recently provided a Shiny app with the same set of functionalities as CellChat Explorer and the advantage of local installation through a docker container.
In the future, we plan to further extend the input options to other published tools. This might, on the one hand, facilitate internal comparison and validation of biological conclusions drawn from interactions predicted by different methods and, on the other hand, give more exibility to the end-user.
In conclusion, InterCellar empowers lab-scientists to interactively analyze cell-cell interactions without programming skills. Moreover, it highlights data-driven patterns of cell-cell communication, by providing several graphical outputs and the possibility to link cellular interactions to biological functions. We believe InterCellar will contribute to accelerate the generation of hypotheses and deepen the understanding of the cellular cross-talk in various pathological and physiological systems.

Implementation
InterCellar is implemented as an R/Shiny application and structured as an R package, which is available on Bioconductor at https://bioconductor.org/packages/InterCellar/. We used the R package golem to ease the development of a robust Shiny application, as well as multiple R packages to build the user interface (shiny, shinydashboard, shinyFiles, shinycssloaders, shinyFeedback, shinyalert, htmltools and htmlwidgets). Data handling is achieved through custom R functions that build upon public R/Bioconductor packages, such as readxl, utils, DT, plyr, dplyr, tidyr and tibble. The functional annotation of InterCellar's function-verse is performed through two R packages that query functional databases, namely biomaRt 30 and graphite 31 . Lastly, graphical outputs are generated by custom R functions based upon R/Bioconductor packages such as ggplot2, plotly, circlize 32 , fmsb, umap, visNetwork, igraph, dendextend 33 , factoextra, colourpicker, scales and grDevices.

Input tools and preprocessing
InterCellar's input data consists of a pre-computed dataset of predicted cell-cell interactions (CCI-data). The application accepts either the output of supported tools (CellPhoneDB v2 and SingleCellSignalR) or of custom analyses performed by the user (for example, as in Kumar et al. 15 ). In the rst case, InterCellar automatically parses the output generated by the supported tool, requiring no further data manipulation by the user. In the second case, the custom input data must contain relevant information structured as a table with the following columns: (i) interaction pairs, containing the names of two molecular components that participate in the interaction (e.g., compA_compB); (ii) communicating cell clusters, containing names or numbers of the cell populations (divided into two columns, e.g., clustA and clustB); (iii) molecular types, either ligand or receptor (in two columns, corresponding to compA and compB); and (iv) a numeric value representing a score for each interaction (e.g., average expression of compA_compB over clustA and clustB).
For both input options, InterCellar will preprocess the data to generate a standardized dataset. These preprocessing steps involve: (1) mapping int-pairs to the associated gene symbols; (2) annotation of the molecular type of each int-pair by combining the provided information with a manually curated set of intpairs (available at https://github.com/martaint/InterCellar-reproducibility); and (3) reordering int-pairs listed as receptor-ligand to ligand-receptor, to ensure consistency in the de nition of the communication ow. To this date, InterCellar supports human genes; the analysis of cell-cell communication obtained from other species can be performed by converting gene names to human orthologs.

Filtering and customization
We implemented multiple ltering options to enable a exible interactive analysis of cell-cell communication. In particular, we distinguish between two types of lters: those applied to the input data or lters applied to visualization options. In the rst case, the user can (i) remove entire clusters from the analysis, for example in the event of unknown or poorly de ned cell populations (in cluster-verse); (ii) re ne int-pairs selecting either a minimum interaction score or a maximum p-value (in cluster-verse); (iii) select int-pairs by annotation source (as in the case of CellPhoneDB) or by speci city (as computed by SingleCellSignalR) (in gene-verse). These three ltering options are applied to the input dataset, therefore the following analyses will be performed on the data subset. The second ltering type concerns graphical outputs and includes features to customize visualizations. These options are only applied to the selected plot and will not affect the underlying dataset. For example, a dot plot visualization available in the geneverse displays int-pair/cluster-pair couplets based on the user's selection of int-pairs of interest. The user can further select a subset of clusters to consider, as well as change the color scale of the dot plot.

Functional annotation of interaction pairs
InterCellar provides interactive annotation of int-pairs with the aim to link biological functions and pathways to cell-cell communication. Speci cally, we use two Bioconductor packages, biomaRt 30 and graphite 31 , to automatically query existing databases of functional terms: Gene Ontology 34,35 , KEGG 19 , Reactome 36 , Biocarta 37 , PID:NCI-Nature 38 , Panther 39 and PharmGKB 40 . InterCellar's functional annotation (available in the function-verse) is required before proceeding to the de nition of int-pair modules. Importantly, the annotation of a certain functional term to an int-pair is ful lled only when all components (e.g., ligand and receptor genes) of the int-pair are enriched by that functional term. Although the user can freely choose which functional sources to consider when performing the annotation, we suggest including as many databases as possible to maximize the number of int-pairs which will be annotated by at least one term. Only these annotated pairs are further considered in the de nition of int-pair modules.

De nition of interaction-pair modules
De ning int-pair modules requires three key, user-driven decisions: (1) choice of the viewpoint cell cluster, representing the cell type of interest; (2) selection of a communication ow, among directed-outgoing (in which the viewpoint cluster sends ligands to other clusters), directed-incoming (in which the viewpoint cluster expresses receptors) and undirected (as in the case of receptor-receptor pairs); and (3) de nition of the number of int-pair modules to consider. In particular, InterCellar subsets the result of the functional annotation (i.e., a binary matrix of int-pairs by functional terms) according to steps (1) and (2). This ltered matrix is given as input to a dimensionality reduction algorithm called UMAP (Uniform Manifold Approximation and Projection 22 ), which calculates a 2D-embedding of the interaction pairs considered (using cosine similarity as metric to compute distances between data points). Thus, the low-dimensional embedding re ects the similarity of int-pairs based on their functional pro les. Finally, using these UMAP coordinates, a hierarchical clustering de nes modules of int-pairs that share similar functional pro les (with euclidean distance and ward.D2 clustering algorithm). The UMAP as well as a further dendrogram result of hierarchical clustering can be visualized and downloaded. In addition, two plots provide guidance to the choice of the optimal number of int-pair modules: (i) total within-cluster sum of squares (WSS) (i.e., elbow method) and (ii) average silhouette width. Both metrics are computed and visualized using functions provided in the R package factoextra. The rst metric indicates the compactness of clusters, by looking at intra-cluster variation. Speci cally, one seeks to minimize the total WSS, represented as a function of the number of clusters. The elbow method, in particular, de nes the optimal number of clusters by looking at a bend (i.e., the elbow) in the plot. We automated this process by using the R package akmedoids, which provides a function to nd the maximum curvature. The average silhouette width evaluates how well each data point lies within its cluster and, on the contrary, should be maximized. By default, InterCellar uses the optimal number of modules computed by the elbow method to generate UMAPs and dendrograms; this value can however be freely changed by the user. When the total number of unique int-pairs considered (by choosing viewpoint and ow) is less or equal to 10, only one module is de ned to prevent noisy results due to insu cient data points.

Signi cant functional terms
To facilitate the biological interpretation of the de ned int-pair modules, InterCellar calculates an empirical p-value for each annotated functional term indicating statistical signi cance of a term to a certain int-pair module. To this end, we implemented a one-sided permutation test which considers as test statistic the proportion of int-pairs annotated to each function, per module (called module ratio). By randomly shu ing the module assignment of each int-pair (1,000 times, without replacement) and calculating the test statistic, we generate the distribution under the null hypothesis. Finally, we compute for each function the proportion of the module ratios which are higher or equal to the actual module ratio, thus obtaining an empirical p-value for the speci city of a function to a certain int-pair module.
Functional terms are then ranked by p-values, and the threshold for signi cance can be chosen by the user (0.05 by default).

COVID-19
This dataset comprises single-cell RNA sequencing data of 24 total patients, divided into 5 healthy controls, 8 COVID-19 moderate cases and 11 COVID-19 critical cases. A total number of 32 samples were derived from nasopharyngeal protected specimen brush and bronchial lavage. Classi cation in moderate and critical cases was performed by the original authors following WHO guidelines. Preprocessed, normalized data were retrieved from Chua et al. 6 , as well as cell type assignment. We removed two cell clusters whose label assignment was poorly de ned, namely "unknown epithelial" and "outlier epithelial", corresponding to ~1.5% of the total number of cells. Moreover, moderate and critical datasets were randomly subsampled to 10,000 cells each (without losing any cell label), while for the control dataset we retained all cells, corresponding to a total of 2,966 cells. Finally, for each of the three datasets, we run CellPhoneDB v2 statistical analysis using default parameters.        Visualization options are divided into three sub-categories to account for speci c features implemented with regards to cell clusters, genes and functional pathways. Automatic download of gures and tables is considered separately.

Supplementary Files
This is a list of supplementary les associated with this preprint. Click to download.