A table of gene counts (matrix of UMI counts for each gene per cell) is required for ICARUS. Datasets with single or multiple samples/replicates may be loaded.

Select between 3 methods of data input:

1. Tab delimited table with cells as columns and gene features as rows (see example A). For datasets with multiple samples, each sample can be denoted by an identifier separated by a underscore (Please do not include underscores for single sample datasets).
2. 10X data files (barcodes.tsv, features.tsv and matrix.mtx, refer to C).
3. Saved Seurat object (RDS file, refer to D)

This step allows the user to load a second dataset that can be integrated with the first dataset. The identical load options are available as previously (refer to Load your data tab). Datasets with single or multiple samples/replicates may be loaded.

Select between 3 methods of data input:

1. Tab delimited table with cells as columns and gene features as rows (see example A). For datasets with multiple samples, each sample can be denoted by an identifier separated by a underscore (Please do not include underscores for single sample datasets).
2. 10X data files (barcodes.tsv, features.tsv and matrix.mtx, refer to C).
3. Saved Seurat object (RDS file, refer to D)

A quality control (QC) step is crucial to remove low quality cells that may arise from cell damage or systematic errors in library preparation (PCR amplification and reverse transcription errors). If not removed, these low-quality cells can distort results in downstream analyses including formation of spurious clusters due to similarities in the damaged-induced expression profiles between otherwise distinct subpopulations.

In this step the user may apply thresholds to commonly used quality control metrics including:

1. Number of unique genes per cell - Low quality cells will have very few genes, cell multiplets will have a high gene count.
2. Number of UMIs per cell - Low quality cells will have very few UMIs, cell multiplets will have a high UMI count.
3. Mitochondrial percentage - Low quality cells will contain high mitochondrial contamination.
4. Ribosomal percentage - Ribosomal percentage will differ depending on the starting cell types, this metric can be useful to removal of cells that do not fall within the expected range.

Cell doublets (multiplets) may arise during scRNA-seq library preparation (i.e., partitioning of multiple cells into a single droplet for droplet-based microfluidics methods). These doublets will exhibit a loss of single cell status and can compromise downstream analysis by creating spurious intermediate populations or transition states.

ICARUS incorporates the DoubletFinder R package to detect and remove doublets. DoubletFinder performs the following:

1. Simulation of "artificial doublets" from existing scRNA-seq data.
2. Merge simulated artificial doublets with real data and perform dimensionality reduction.
3. A k-nearest neighbour graph is developed, and each cell is scored based on its proximity to artificial doublets. The highest scoring cells are assigned as real doublets.
4. The user has the option to visualise and remove these doublets.

Please note doublet removal may take a long time depending on the size of the input dataset(s). Please save the output and use the load function for repeat analysis.

A quality control (QC) step is crucial to remove low quality cells that may arise from cell damage or systematic errors in library preparation (PCR amplification and reverse transcription errors). If not removed, these low-quality cells can distort results in downstream analyses including formation of spurious clusters due to similarities in the damaged-induced expression profiles between otherwise distinct subpopulations.

In this step the user may apply thresholds to commonly used quality control metrics including:

1. Number of unique genes per cell - Low quality cells will have very few genes, cell multiplets will have a high gene count.
2. Number of UMIs per cell - Low quality cells will have very few UMIs, cell multiplets will have a high UMI count.
3. Mitochondrial percentage - Low quality cells will contain high mitochondrial contamination.
4. Ribosomal percentage - Ribosomal percentage will differ depending on the starting cell types, this metric can be useful to removal of cells that do not fall within the expected range.

Cell doublets (multiplets) may arise during scRNA-seq library preparation (i.e., partitioning of multiple cells into a single droplet for droplet-based microfluidics methods). These doublets will exhibit a loss of single cell status and can compromise downstream analysis by creating spurious intermediate populations or transition states.

ICARUS incorporates the DoubletFinder R package to detect and remove doublets. DoubletFinder performs the following:

1. Simulation of "artificial doublets" from existing scRNA-seq data.
2. Merge simulated artificial doublets with real data and perform dimensionality reduction.
3. A k-nearest neighbour graph is developed, and each cell is scored based on its proximity to artificial doublets. The highest scoring cells are assigned as real doublets.
4. The user has the option to visualise and remove these doublets.

Please note doublet removal may take a long time depending on the size of the input dataset(s). Please save the output and use the load function for repeat analysis.

In single cell RNA-seq analyses, individual cells may be characterised by its transcriptomic expression profile comprised of thousands of genes. If the expression level of each gene represents a dimension of data, each cell can be allocated a location in the high-dimensional expression space. The purpose of dimensionality reduction is to enable the analysis and visualisation of high-dimensional data in a low dimensional plane whilst retaining the major properties of the original data. This is possible as biological processes encompass multiple genes (i.e., genes within a pathway) and therefore, multiple features can be compressed into a single dimension.

The following are performed at this step:

1. Either log normalization of the data by a scale factor of 10,000 or normalization by SCTransform. A Normalization step ensures technical biases or noise present from sequencing are kept to a minimum.
2. A linear transformation is applied to scale expression levels so that the mean expression of genes across cells equates to zero and the variance across cells equates to one. This step gives equal weight in downstream analyses, so that highly-expressed genes do not dominate.
3. Dimensionality reduction using Principal Component Analysis (PCA) with the user defined number of variable features to incorporate. Seurat performs dimensionality reduction on a set of highly variable genes to improve computational speed and highlight biological signals in single-cell datasets, this number may be altered in the NEW box below.
4. An option to impute dropouts (false zeros in the dataset due to low amounts of mRNA in individual cells resulting in insufficient mRNA capture) is available. The Adaptively-thresholded low rank approximation (ALRA) method to impute dropouts is used.

In single cell RNA-seq analyses, individual cells may be characterised by its transcriptomic expression profile comprised of thousands of genes. If the expression level of each gene represents a dimension of data, each cell can be allocated a location in the high-dimensional expression space. The purpose of dimensionality reduction is to enable the analysis and visualisation of high-dimensional data in a low dimensional plane whilst retaining the major properties of the original data. This is possible as biological processes encompass multiple genes (i.e., genes within a pathway) and therefore, multiple features can be compressed into a single dimension.

The following are performed at this step:

1. Either log normalization of the data by a scale factor of 10,000 or normalization by SCTransform. A Normalization step ensures technical biases or noise present from sequencing are kept to a minimum.
2. A linear transformation is applied to scale expression levels so that the mean expression of genes across cells equates to zero and the variance across cells equates to one. This step gives equal weight in downstream analyses, so that highly-expressed genes do not dominate.
3. Dimensionality reduction is performed using Principal Component Analysis (PCA) with the user defined number of variable features. Seurat performs dimensionality reduction on a set of highly variable genes to improve computational speed and highlight biological signals in single-cell datasets, this number may be altered in the NEW box below.
4. Data integration with either Canonical Correlation Analysis (CCA), Reciprocal PCA (RPCA) or Harmony methodologies. For CCA and RPCA methods, the number of k-anchors (strength of integration) can be adjusted, please increase the k-anchors parameter (default value = 5) for samples where integration of certain cell types are not aligned.

The outcome of clustering is to group together cells with similar expression profiles for human interpretation.

ICARUS employs Seurat's graph-based community detection clustering approach. In brief, a graph of k-nearest neighbours (k-NN) is formed between cells in high dimensional space where each cell is a node that is connected to its nearest cells. The edges of each connection are then weighted based on its similarities to neighbouring cells and a community detection algorithm is utilised to define clusters.

Clustering parameters including number of dimensions, number of k-nearest neighbours for graph construction and choice of clustering algorithm (either Louvain, SLM or Leiden) may be adjusted in this step. UMAP parameters (nearest neighbours and minimum distance) and t-SNE parameters (perplexity and number of iterations) can also be adjusted.

Please note clustering may take a long time depending on the size of the input dataset(s). Please save the output and use the load function for repeat analysis.

Data correction aims to remove further technical and biological confounders. For example, cell clusters may be formed largely due to cell cycle state rather than differences in underlying biology. The choice whether to correct for these effects will depend on the intention of the downstream analysis.

This step allows the regression of the following:

1. Cell cycle genes: Expression of cell cycle phase canonical markers (taken from Triosh et al., 2016).
2. User gene(s) of interest: Average normalised gene count for each selected gene are computed and regressed from the data (using the vars.to.regress function in Seurat::ScaleData).
3. User gene pathway of interest: Average normalised gene count for all genes in the gene pathway (GSEA/Msigdb database) are computed for each cell and regressed from the data (using the vars.to.regress function in Seurat::ScaleData).

Please note data correction may take a long time depending on the size of the input dataset(s). Please save the output and use the load function for repeat analysis.

Marker genes for each cluster may be analysed and compared against known cell markers to assign a cluster cell type. ICARUS enables users to add custom cell cluster labels or annotate clusters against either the sctype reference set or the singleR reference set.

The sctype reference dataset compiled cell markers from CellMarker and PanglaoDB databases. These include markers for the brain, pancreas, immune system, liver, adrenal, heart, eye, kidney, lung, muscle, placenta, spleen, stomach and thymus.

The singleR R package incorporates cell markers from the following databases:
Brain

1. Darmanis Brain Data: scRNA-seq of 8 adult and 4 fetal human brain temporal lobe tissue (Darmanis et al., 2015).
2. Zhong Prefrontal Cortex Data: scRNA-seq of developing human embryonic prefrontal cortex from gestational weeks 8 to 26 (Zhong et al., 2018).

Immune system

1. Blueprint Encode: scRNA-seq of 259 bulk RNA-seq samples generated by Blueprint and ENCODE from pure populations of stroma and immune cells (Martens and Stunnenberg, 2013; The ENCODE Consortium, 2012).
2. Database of Immune Cell Expression (DICE): scRNA-seq of 1561 bulk RNA-seq samples generated by DICE from pure populations of human immune cells. (Schmiedel et al., 2018).
3. Immunologic Genome Project (ImmGen): scRNA-seq of 830 microarray samples generated by ImmGen from pure populations of murine immune cells (http://www.immgen.org/).
4. Monaco Immune Data: scRNA-seq of 114 bulk RNA-seq samples of sorted immune cell populations that can be found in GSE107011 (Monaco et al., 2019).
5. Novershtern Hematopoietic Data: scRNA-seq of 211 bulk human microarray samples of sorted hematopoietic cell populations that can be found in GSE24759 (Novershtern et al., 2011).

Pancreas

1. Baron Pancreas Data: scRNA-seq of pancreatic cells from 4 human donors and 2 mouse strains (Baron et al., 2016).
2. Lawlor Pancreas Data: scRNA-seq of 5 non-diabetic and 8 type 2 diabetic human islet samples (Lawlor et al., 2017).
3. Muraro Pancreas Data: scRNA-seq of pancreatic cells from 4 human donors (Muraro et al., 2016).
4. Segerstolpe Pancreas Data: scRNA-seq of human islet cells from 6 healthy and 4 type 2 diabetic donors (Segerstolpe et al., 2016).

Multiple organs

1. Mouse RNA-seq: 358 bulk RNA-seq samples of sorted cell populations (Benayoun et al., 2019).
2. Human Primary Cell Atlas: scRNA-seq of 713 microarray samples from the Human Primary Cell Atlas (Mabbott et al., 2013).
3. He Organ Atlas: scRNA-seq of 15 tissue organs of one adult donor. Tissue organs include bladder, blood, common bile duct, oesophagus, heart, liver, lymph node, bone marrow, muscle, rectum, skin, small intestine, spleen, stomach and trachea (He et al., 2020).

Please note cluster labelling using reference marker databases may take a long time depending on the size of the input dataset(s). Please save the output and use the load function for repeat analysis.

This step allows the user to visualize normalised counts for gene(s) or gene pathway (GSEA/Msigdb gene pathway including KEGG, Reactome, WikiPathways and Gene Ontology) across cell clusters.

Many biological networks contain patterns of gene expression that are evolutionary conserved. These patterns can be characterised and grouped into hierarchical co-expression modules which can infer biological functionality. A method of detecting these gene Co-Expression modules is MEGENA (Multiscale Embedded Gene Co-Expression Network Analysis) which employs a Planar Maximally filtered graph to extract significant gene interactions and networks of co-expression modules.

MEGENA identifies gene co-expression networks through:

1. Identification of significant gene interactions (Planar filtered network (PFN) construction through embedding on a topological sphere).
2. Multiscale clustering analysis to group similar and dissimilar gene modules.
3. Multiscale hub analysis to detect connected hubs of individual clusters.

ICARUS performs MEGENA on a set of Seurat highly variable genes (see Dimensionality Reduction tab) or differentially expressed genes to generate a set of hierarchically ordered co-expression modules, where larger modules progressively branch into smaller submodules.

A cell’s transcriptional state may be characterised by gene regulatory networks (GRNs) that are formed from transcription factors and cofactors that regulate each other and their downstream gene targets. ICARUS employs the SCENIC R package to characterise cell cluster/cell type specific GRNs using a set of variable genes (seurat variable features, please refer to the Dimensionality Reduction tab) or user computed differentially expressed genes (refer to the differential expression and custom differential expression tabs).

The SCENIC workflow consists of:

1. Identification of coexpression modules with GENIE3.
2. Cis-regulatory transcription factor binding motif analysis of co-expressed transcription factors and genes. Motifs in the promoter of genes (up to 500bp upstream the TSS) and in the 20kb around the TSS (+/- 10kb) are scored. Co-expression modules with significant motif enrichment are retained and indirect targets are removed. These modules are termed regulons.
3. The activity of these transcription factor regulated gene modules (regulons) are scored across cell clusters/cell types using the AUCell algorithm.

WARNING: ICARUS by default limits regulon construction to 100 cells per cluster/cell type to improve computational speed. The user may choose to increase the number of cells analysed. However, these analyses could take hour to days depending on the size of the analysis.

Genome wide association studies have identified several loci in various genes that are associated with a trait of interest. The expression profiles of cell clusters can be compared against these GWAS loci to identify potential causal cell types underlying complex traits.

ICARUS employs the multi-marker analysis of genomic annotation (MAGMA) methodology to identify increased linear association between cluster derived gene sets and GWAS traits. The user may either upload their own GWAS summary statistics or select biological traits from hundreds of GWAS datasets curated by Neurogenomics/MAGMA_Files_Public which includes gene level association statistics for 10b upstream and 1.5kb downstream of each gene.

The expression profiles of cells change during development, in response to stimuli and throughout life. Trajectory analysis aims to determine the sequence of gene expression changes in your dataset. ICARUS employs the Monocle3 algorithm to graph cells according to their progress in pseudotime (a measure of how much progress an individual cell has made, i.e., cell differentiation). Pseudotime estimates the amount of transcriptional change a cell undergoes from its beginning to end states.

Please use the lasso tool that can be found at the top right hand corner of the interactive plot to select your groups of cells.

Cells interact and communicate with each other. The coordination of cell-cell communication is vital for many biological processes including cell cycle activation, initiation of apoptosis and activation of cellular development and differentiation.

ICARUS employs the CellChat tool to quantitatively infer and analyse intercellular communication networks from single cell RNA-seq data. CellChat incorporates a manually curated comprehensive database of ligand-receptor pairs, soluble agonists/antagonists and stimulatory/inhibitory membrane bound co-receptors to infer cell-cell communication interactions based on social network analysis tools, pattern recognition methods and manifold learning approaches.

The CellChat workflow consists of:

1. Identification of differentially expressed ligands and receptors genes for each cell cluster/cell type.
2. A communication probability (interaction strength) is computed by modelling ligand-receptor interactions using the law of mass action on average expression values of a ligand in one cell cluster/cell type and a receptor of another cell cluster/cell group.
3. Significant interactions are identified using a permutation test that randomly permutes cell group labels and recomputes the interaction probability.

Multimodal analysis refers to the simultaneous measurements of several data types from the same cell. For example, methodologies including CITE-seq, cell hasing oligos and 10X multiome kit allow measurements of single cell transcriptomes and cell-surface proteins of the same cell. ICARUS enables a side by side visualisation of multi-modal data to enable easy comparisons and analysis.

ICARUS requires a matrix of multimodal data gene counts. This may be uploaded as a:

1. Tab delimited table with cells as columns and gene features as rows (see example A). Please make sure the cell names (column names) match the loaded scRNA-seq cell names.
   
2. 10X data files (barcodes.tsv, features.tsv and matrix.mtx, refer to B).
3. Saved Seurat object (RDS file, refer to C).

ICARUS enables users to perform pairwise differential expression tests between clusters of interest or between samples within cell cluster(s) and output a list of differentially expressed genes. The list of differentially expressed genes are then probed for enriched gene pathways (gene set enrichment analysis) to reveal potential affected biological pathways. Gene set enrichment and visualisation performed using ClusterProfiler and ReactomePA R packages.

Several outputs are produced in this step including:

1. Table of differentially expressed genes.
2. Table of enriched pathways (Gene Ontology/WikiPathways/Reactome/KEGG pathways).
3. Volcano plot of differentially expressed genes.
4. Visualisation of log2 fold changes for specific gene pathway.

This step allows the user to investigate genes that are differentially expressed between user selected cells. Pairwise differential expression tests are computed between selected cell groups and a list of differentially expressed genes is produced. Gene set enrichment analysis of significantly differentially expressed genes may reveal affected gene pathways. Gene set enrichment and visualisation performed using ClusterProfiler and ReactomePA R packages.

Several outputs are produced in this step including:

1. Table of differentially expressed genes.
2. Table of enriched pathways (Gene Ontology/WikiPathways/Reactome/KEGG pathways).
3. Volcano plot of differentially expressed genes.
4. Visualisation of log2 fold changes for specific gene pathway.

Please use the lasso tool that can be found at the top right hand corner of the interactive plot to select your groups of cells. Press COMPARE to compare the two groups of lasso selected cells.

Additional graphical outputs for gene set enrichment analysis. The outputs of this tab use the list of differentially expressed genes computed in the "Differential Expression" tab.

Visualisations of enriched terms in the form of:

1. Dot plot: Enriched terms are plotted in order of gene ratio (number of enriched genes in gene pathway/total number of genes in gene pathway).
2. Gene Concept Network: Visualisation of complex association and interactions between genes and enriched terms as a network.
3. Enrichment map: Network of enriched terms with edges connecting overlapping gene pathways. Mutually overlapping gene pathways are clustered together allowing easy identification of related pathways.
4. View specific pathways: View selected reactome pathway.

The Drug Gene Interaction Database (DGIdb) curates information on druggable genes from publications, known affected pathways (Gene Ontology, Human Protein Atlas, IDG) and publicly available databases (DrugBank, PharmGKB, Chembl, Drug Target Commons and TTD). ICARUS offers the user the option to query differentially expressed genes (refer to differential expression tab) against DGIdb v4.0 and returns a list of potential drug targets with information of interaction type (i.e., inhibitor, agonist, blocker), interaction claim source (i.e., PharmGKB, ChemblInteractions, etc), interaction group score (score takes into account number of drug and gene partners and number of supporting publications) and the relevant pubmed reference if available.

Identification of drug-gene targets can facilitate targeted perturbations of key molecular pathways. In the context of dysregulation or disease, this could provide an avenue towards a therapeutic opportunity through repurposed drugs.

Additional graphical outputs for gene set enrichment analysis. The outputs of this tab use the set of differentially expressed genes computed in the "Custom Differential Expression" tab.

Visualisations of enriched terms in the form of:

1. Dot plot: Enriched terms are plotted in order of gene ratio (number of enriched genes in gene pathway/total number of genes in gene pathway).
2. Gene Concept Network: Visualisation of complex association and interactions between genes and enriched terms as a network.
3. Enrichment map: Network of enriched terms with edges connecting overlapping gene pathways. Mutually overlapping gene pathways are clustered together allowing easy identification of related pathways.
4. View specific pathways: View selected reactome pathway.

The Drug Gene Interaction Database (DGIdb) curates information on druggable genes from publications, known affected pathways (Gene Ontology, Human Protein Atlas, IDG) and publicly available databases (DrugBank, PharmGKB, Chembl, Drug Target Commons and TTD). ICARUS offers the user the option to query differentially expressed genes (refer to custom differential expression tab) against DGIdb v4.0 and returns a list of potential drug targets with information of interaction type (i.e., inhibitor, agonist, blocker), interaction claim source (i.e., PharmGKB, ChemblInteractions, etc), interaction group score (score takes into account number of drug and gene partners and number of supporting publications) and the relevant pubmed reference if available.

Identification of drug-gene targets can facilitate targeted perturbations of key molecular pathways. In the context of dysregulation or disease, this could provide an avenue towards a therapeutic opportunity through repurposed drugs.