Biología de Sistemas Unizar: Sox2 Enrichment
Andrés Gordo Ortiz
March 15, 2024
1 What can I find in this Report?
This analysis performs a fully crafted Transcriptomics Protocol, from read mapping, filtering and normalisation, up to Differential Gene Expression or Ontology enrichment
1.1 Data Integration
- Following count matrix import from Arcs4, the TxImport package was utilized to import Kallisto outputs into the R environment.
1.2 Data Summarization
- Annotation data from Biomart was employed to summarize the data from transcript-level to gene-level, providing a comprehensive view of the gene expression landscape.
2 Preprocessing
2.1 Filtering & normalisation
Filtering was carried out to remove lowly expressed genes. Genes with less than 1 count per million (CPM) in at least 615 or more samples were filtered out. This is done to make sure that the low expression is due to a lack of significance across all conditions, and not an intrinsic property of one of them. This filtering reduced the number of genes from 67186 to 12988. Normalisation of samples was performed with edgeR, using Trimmed Mean of M-value or TMM. This method is based on the assumption that most genes are not differentially expressed. It calculates a scaling factor for each sample, which is the median of the ratio of each gene’s expression to the geometric mean of all samples. This scaling factor is then used to normalise the expression of each gene in each sample. The normalised expression is then expressed as log2(CPM). The normalisation step is important because it allows for the comparison of expression between samples.
According to fastQC results, a total of 0 samples were removed due to low quality. This is NaN% of the total reads. Filtered out samples are: .
2.2 Filtered and Normalised data
Scatter plot showing the divergence in the average expression of every gene between the two conditions
The Table 1 includes expression data for 12988 genes. You can sort and search the data directly from the table.
3 Principal Component Analysis
3.1 PCA
Principal Component Analysis plot coloured by condition. The 95% confidence interval is shown if enough samples are analysed.
4 Volcano plot
An important note: all the contrasts, design matrix and differential expression analysis are made according to this formula (SOX UP - SOX DOWN). This means that a logfold change > 0 shows that the first condition is upregulated, whereas a logfold change < 0 shows downregulation.
Volcano plot of differencial gene expression.The lines delimit the log-fold change and adjusted p-value (FDR) thresholds.
5 MA Plot
The plot is an MA plot colored by False Discovery Rate (FDR). It shows the relationship between average expression level (AveExpr) and log fold change (logFC) for genes or features. Points are colored blue for FDR passed, or red for FDR not passed
6 Table of DEGs
To identify differentially expressed genes, precision weights were first applied to each gene based on its mean-variance relationship using VOOM, then data was normalized using the TMM method in EdgeR. Linear modeling and bayesian stats were employed via Limma to find genes that were up- or down-regulated by 3-logfold or more, with a false-discovery rate (FDR) of 10^{-100}.
7 Heatmaps and Modules
Pearson correlation was used to cluster 4529 differentially expressed genes, which were then represented as heatmap with the data scaled by Zscore for each row. On the other hand, Spearman correlation was used to cluster all samples. Two modules can be seen, being either up- or down-regulated depending on the condition of the samples.
7.1 Complete
7.2 Upregulated in SOX UP
7.3 Downregulated in SOX UP
8 GO enrichment
Gene Ontology enrichment for the 12988 genes differentially expressed.
8.1 Upregulated genes in SOX UP
8.1.1 Manhattan Plot
Manhattan plot of GO enrichment for upregulated genes
8.1.2 Table
8.2 Downregulated genes in SOX UP
8.2.1 Manhattan Plot
Manhattan plot of GO enrichment for downregulated genes
8.2.2 Table
9 GSEA
The 34550 gene sets in the Human Molecular Signatures Database (MSigDB) are divided into 9 major collections, and several subcollections. In this analysis, we will use the C2 (Curated Genes) and H (Hallmark genes) collections, which contains gene sets that represent canonical pathways, gene ontology, and other gene sets derived from knowledge in the literature. The top 15 genes set for each conditions and collection will have their own Gene Set Enrichment Plot. The table will contain all the results for each gene set. A sample plot is shown below.
9.1 C2 Collection
9.1.1 GSEA Plot
Sample of a Gene Set Enrichment Plot from the C2 Collection
9.1.2 Table
9.2 Hallmark Collection
9.2.1 GSEA Plot
Sample of a Gene Set Enrichment Plot from the Hallmark Collection
9.2.2 Table
9.3 Bubble Plot
Finally, a Bubble Plot is produced for the top 20 gene sets. The size of the bubble is proportional to the number of genes in the gene set, the color is proportional to the NES, and the alpha is proportional to the -log10(p.adjust). The plot is shown below. NES means normalized enrichment score, and it is the primary statistic for ranking genes in a GSEA analysis. It represents the degree to which a gene set is overrepresented at the top or bottom of a ranked list of genes. The score is normalized to account for differences in gene set size and in correlations between gene set members and the expression dataset. The NES score can be either >0, which means is enriched in the SOX UP condition, or <0, which means it is enriched in the SOX DOWN condition. The table provided shows the absolute score beacause it has already been classified in one the two conditions.
9.3.1 C2 Collection
Bubble plot of the top 20 gene sets within the C2 Collection
9.3.2 Hallmarks Collection
Bubble plot of the top 20 gene sets within the H Collection
10 Session info
The output from running ‘sessionInfo’ is shown below and details all packages and version necessary to reproduce the results in this report.
## R version 4.3.1 (2023-06-16 ucrt)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 11 x64 (build 22621)
##
## Matrix products: default
##
##
## locale:
## [1] LC_COLLATE=English_United Kingdom.utf8
## [2] LC_CTYPE=English_United Kingdom.utf8
## [3] LC_MONETARY=English_United Kingdom.utf8
## [4] LC_NUMERIC=C
## [5] LC_TIME=English_United Kingdom.utf8
##
## time zone: Europe/Madrid
## tzcode source: internal
##
## attached base packages:
## [1] stats4 stats graphics grDevices datasets utils methods
## [8] base
##
## other attached packages:
## [1] IsoformSwitchAnalyzeR_2.2.0 pfamAnalyzeR_1.2.0
## [3] sva_3.50.0 genefilter_1.84.0
## [5] mgcv_1.9-1 nlme_3.1-164
## [7] satuRn_1.10.0 DEXSeq_1.48.0
## [9] DESeq2_1.42.0 SummarizedExperiment_1.32.0
## [11] MatrixGenerics_1.14.0 BiocParallel_1.36.0
## [13] enrichplot_1.22.0 msigdbr_7.5.1
## [15] clusterProfiler_4.10.0 gprofiler2_0.2.2
## [17] GSVA_1.50.0 GSEABase_1.64.0
## [19] graph_1.80.0 annotate_1.80.0
## [21] XML_3.99-0.16.1 heatmaply_1.5.0
## [23] viridis_0.6.5 viridisLite_0.4.2
## [25] RColorBrewer_1.1-3 gplots_3.1.3.1
## [27] plotly_4.10.4 gt_0.10.1
## [29] DT_0.31 svglite_2.1.3
## [31] cowplot_1.1.3 matrixStats_1.2.0
## [33] edgeR_4.0.16 limma_3.58.1
## [35] EnsDb.Hsapiens.v86_2.99.0 ensembldb_2.26.0
## [37] AnnotationFilter_1.26.0 GenomicFeatures_1.54.3
## [39] AnnotationDbi_1.64.1 Biobase_2.62.0
## [41] GenomicRanges_1.54.1 GenomeInfoDb_1.38.6
## [43] IRanges_2.36.0 S4Vectors_0.40.2
## [45] BiocGenerics_0.48.1 tximport_1.30.0
## [47] lubridate_1.9.3 forcats_1.0.0
## [49] stringr_1.5.1 dplyr_1.1.4
## [51] purrr_1.0.2 readr_2.1.5
## [53] tidyr_1.3.1 tibble_3.2.1
## [55] ggplot2_3.4.4 tidyverse_2.0.0
## [57] beepr_1.3 knitr_1.45
## [59] readxl_1.4.3 openxlsx_4.2.5.2
## [61] tinytex_0.49 rmarkdown_2.25
##
## loaded via a namespace (and not attached):
## [1] progress_1.2.3 locfdr_1.1-8
## [3] Biostrings_2.70.2 HDF5Array_1.30.1
## [5] vctrs_0.6.5 digest_0.6.34
## [7] png_0.1-8 registry_0.5-1
## [9] ggrepel_0.9.5 renv_1.0.3
## [11] MASS_7.3-60.0.1 reshape2_1.4.4
## [13] foreach_1.5.2 httpuv_1.6.14
## [15] qvalue_2.34.0 withr_3.0.0
## [17] xfun_0.42 ggfun_0.1.4
## [19] ellipsis_0.3.2 survival_3.5-8
## [21] memoise_2.0.1 tximeta_1.20.3
## [23] gson_0.1.0 systemfonts_1.0.5
## [25] ragg_1.2.7 tidytree_0.4.6
## [27] gtools_3.9.5 pbapply_1.7-2
## [29] prettyunits_1.2.0 KEGGREST_1.42.0
## [31] promises_1.2.1 httr_1.4.7
## [33] restfulr_0.0.15 rhdf5filters_1.14.1
## [35] rhdf5_2.46.1 rstudioapi_0.15.0
## [37] generics_0.1.3 DOSE_3.28.2
## [39] babelgene_22.9 curl_5.2.0
## [41] zlibbioc_1.48.0 ScaledMatrix_1.10.0
## [43] ggraph_2.1.0 polyclip_1.10-6
## [45] ca_0.71.1 GenomeInfoDbData_1.2.11
## [47] SparseArray_1.2.4 interactiveDisplayBase_1.40.0
## [49] xtable_1.8-4 evaluate_0.23
## [51] S4Arrays_1.2.0 BiocFileCache_2.10.1
## [53] hms_1.1.3 irlba_2.3.5.1
## [55] colorspace_2.1-0 filelock_1.0.3
## [57] VennDiagram_1.7.3 magrittr_2.0.3
## [59] later_1.3.2 ggtree_3.10.0
## [61] lattice_0.22-5 shadowtext_0.1.3
## [63] pillar_1.9.0 iterators_1.0.14
## [65] caTools_1.18.2 compiler_4.3.1
## [67] beachmat_2.18.1 stringi_1.8.3
## [69] TSP_1.2-4 dendextend_1.17.1
## [71] GenomicAlignments_1.38.2 plyr_1.8.9
## [73] crayon_1.5.2 abind_1.4-5
## [75] BiocIO_1.12.0 gridGraphics_0.5-1
## [77] locfit_1.5-9.8 graphlayouts_1.1.0
## [79] bit_4.0.5 fastmatch_1.1-4
## [81] textshaping_0.3.7 codetools_0.2-19
## [83] BiocSingular_1.18.0 crosstalk_1.2.1
## [85] bslib_0.6.1 mime_0.12
## [87] splines_4.3.1 Rcpp_1.0.12
## [89] dbplyr_2.4.0 sparseMatrixStats_1.14.0
## [91] HDO.db_0.99.1 cellranger_1.1.0
## [93] blob_1.2.4 utf8_1.2.4
## [95] BiocVersion_3.18.1 fs_1.6.3
## [97] DelayedMatrixStats_1.24.0 ggplotify_0.1.2
## [99] Matrix_1.6-5 statmod_1.5.0
## [101] tzdb_0.4.0 tweenr_2.0.2
## [103] pkgconfig_2.0.3 tools_4.3.1
## [105] cachem_1.0.8 RSQLite_2.3.5
## [107] DBI_1.2.2 fastmap_1.1.1
## [109] scales_1.3.0 grid_4.3.1
## [111] audio_0.1-11 geneplotter_1.80.0
## [113] Rsamtools_2.18.0 AnnotationHub_3.10.0
## [115] sass_0.4.8 patchwork_1.2.0
## [117] BiocManager_1.30.22 snow_0.4-4
## [119] farver_2.1.1 tidygraph_1.3.1
## [121] scatterpie_0.2.1 yaml_2.3.8
## [123] rtracklayer_1.62.0 cli_3.6.2
## [125] webshot_0.5.5 lifecycle_1.0.4
## [127] lambda.r_1.2.4 timechange_0.3.0
## [129] gtable_0.3.4 rjson_0.2.21
## [131] parallel_4.3.1 ape_5.7-1
## [133] jsonlite_1.8.8 seriation_1.5.4
## [135] bitops_1.0-7 bit64_4.0.5
## [137] assertthat_0.2.1 yulab.utils_0.1.4
## [139] zip_2.3.1 futile.options_1.0.1
## [141] highr_0.10 jquerylib_0.1.4
## [143] GOSemSim_2.28.1 lazyeval_0.2.2
## [145] shiny_1.8.0 htmltools_0.5.7
## [147] GO.db_3.18.0 rappdirs_0.3.3
## [149] formatR_1.14 glue_1.7.0
## [151] XVector_0.42.0 RCurl_1.98-1.14
## [153] treeio_1.27.0.002 BSgenome_1.70.2
## [155] futile.logger_1.4.3 gridExtra_2.3
## [157] boot_1.3-28.1 igraph_2.0.2
## [159] R6_2.5.1 SingleCellExperiment_1.24.0
## [161] labeling_0.4.3 Rhdf5lib_1.24.2
## [163] aplot_0.2.2 DelayedArray_0.28.0
## [165] tidyselect_1.2.0 ProtGenerics_1.34.0
## [167] ggforce_0.4.1 xml2_1.3.6
## [169] rsvd_1.0.5 munsell_0.5.0
## [171] KernSmooth_2.23-22 data.table_1.15.0
## [173] htmlwidgets_1.6.4 fgsea_1.28.0
## [175] hwriter_1.3.2.1 biomaRt_2.58.2
## [177] rlang_1.1.3 fansi_1.0.6