Explore SingleCellExperiment using tidyverse

Summary

• retina dataset
  • Explore using cellxgene GUI
  • Download Seurat .rds object
  • Published as Voigt et al., 2019
• Student references
  • https://swirlstats.com
  • https://r4ds.had.co.nz
  • http://bioconductor.org/books/OSCA

Load packages

library( "Seurat" )
library( "scater" )

library( "tidyverse" )

System has not been booted with systemd as init system (PID 1). Can't operate.
Failed to create bus connection: Host is down
running command 'timedatectl' had status 1

Load data

sce <- readRDS( "/home/idies/workspace/c_moor_data/retina/retina.rds" )
sce <- as.SingleCellExperiment( sce )
sce

class: SingleCellExperiment 
dim: 21599 8217 
metadata(0):
assays(1): logcounts
rownames(21599): FAM138A RP11-34P13.7 ... AC145205.1
  AC145205.1-1
rowData names(1): hgnc_gene_symbol
colnames(8217): AAACCCAGTCGGCACT-1 AAACGAACAGGGACTA-1 ...
  TTTGGAGCAGTTTGGT-6 TTTGGAGTCGTTCTAT-6
colData names(18): final_cluster_labels libraryID ...
  cell_type ident
reducedDimNames(1): UMAP
altExpNames(0):

Explore data

UMAP

plotUMAP( sce, colour_by="cell_type" )

Table

table( sce$cell_type )

                amacrine cell                bipolar neuron 
                           98                          2391 
             endothelial cell     foveal cone photoreceptor 
                          254                           236 
                   glial cell               microglial cell 
                         2413                           154 
                pericyte cell peripheral cone photoreceptor 
                          111                            68 
       retina horizontal cell         retinal ganglion cell 
                           82                           326 
             retinal rod cell                       unknown 
                         1590                           494 

Expression

plotExpression( sce, "RHO", "cell_type" ) +
  theme( axis.text.x = element_text( angle=90 ) )

Prepare data.frame

Markers

df_markers <- tribble(
  ~marker, ~cell_type,
  #------/------------
  "RHO", "rod",
  "PDE6A", "rod",
  "PDE6H", "cone",
  "ARR3", "cone",
  "TRPM1", "bipolar",
  "THY1", "RGC",
  "NRN1", "RGC",
  "ONECUT2", "horizontal",
  "NRXN2", "amacrine"
)
goi <- df_markers$marker

Counts

df_counts <- assay( sce, "logcounts" )[goi,] %>%
  t() %>%
  data.frame() %>%
  rownames_to_column( "barcode" )
# df_counts

colData

coi <- c( "cell_type", "celltype", "region", "disease", "sex", "libraryID" )

df_coldata <- colData( sce )[,coi] %>%
  data.frame() %>%
  rownames_to_column( "barcode" )

df_umap <- reducedDim( sce, "UMAP" ) %>%
  data.frame() %>%
  rownames_to_column( "barcode" )

df

df <- df_coldata %>%
  full_join( df_umap, by="barcode" ) %>%
  full_join (df_counts, by="barcode" )
# df

colnames( df )

 [1] "barcode"   "cell_type" "celltype"  "region"    "disease"  
 [6] "sex"       "libraryID" "UMAP_1"    "UMAP_2"    "RHO"      
[11] "PDE6A"     "PDE6H"     "ARR3"      "TRPM1"     "THY1"     
[16] "NRN1"      "ONECUT2"   "NRXN2"    

Tidyverse

UMAP cell_type

ggplot( df ) +
  geom_point( aes( UMAP_1, UMAP_2, color=cell_type ) )

UMAP region

ggplot( df ) +
  geom_point( aes( UMAP_1, UMAP_2, color=region ) )

UMAP RHO

ggplot( df ) +
  geom_point( aes( UMAP_1, UMAP_2, color=RHO ) ) 

Expression RHO

ggplot( df ) +
  geom_violin( aes( cell_type, RHO ), scale="width" )

Expression Facet

df %>%
  tidyr::pivot_longer( goi, "gene" ) %>%
  ggplot() +
  geom_violin( aes( gene, value ), scale="width" ) +
  theme( axis.text.x = element_text( angle=90 ) ) +
  facet_wrap( ~cell_type, scales="free_y" )

Note: Using an external vector in selections is ambiguous.
ℹ Use `all_of(goi)` instead of `goi` to silence this message.
ℹ See <https://tidyselect.r-lib.org/reference/faq-external-vector.html>.
This message is displayed once per session.

Composition

ggplot( df ) +
  geom_bar( aes( libraryID, fill=cell_type ), position="fill" ) +
  theme( axis.text.x = element_text( angle=90 ) )

Table

df %>%
  group_by( cell_type ) %>%
  count()

Two Categories?

ggplot( df ) +
  geom_count( aes( libraryID, cell_type ) ) +
  theme( axis.text.x = element_text( angle=90 ) )

Document analysis

sessionInfo()

R version 4.0.3 (2020-10-10)
Platform: x86_64-pc-linux-gnu (64-bit)
Running under: CentOS Linux 8

Matrix products: default
BLAS:   /home/idies/R/lib64/R/lib/libRblas.so
LAPACK: /home/idies/R/lib64/R/lib/libRlapack.so

locale:
 [1] LC_CTYPE=en_US.UTF-8       LC_NUMERIC=C              
 [3] LC_TIME=en_US.UTF-8        LC_COLLATE=en_US.UTF-8    
 [5] LC_MONETARY=en_US.UTF-8    LC_MESSAGES=en_US.UTF-8   
 [7] LC_PAPER=en_US.UTF-8       LC_NAME=C                 
 [9] LC_ADDRESS=C               LC_TELEPHONE=C            
[11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C       

attached base packages:
[1] parallel  stats4    stats     graphics  grDevices utils    
[7] datasets  methods   base     

other attached packages:
 [1] forcats_0.5.1               stringr_1.4.0              
 [3] dplyr_1.0.5                 purrr_0.3.4                
 [5] readr_1.4.0                 tidyr_1.1.3                
 [7] tibble_3.1.1                tidyverse_1.3.1            
 [9] scater_1.18.6               ggplot2_3.3.3              
[11] SingleCellExperiment_1.12.0 SummarizedExperiment_1.20.0
[13] Biobase_2.50.0              GenomicRanges_1.42.0       
[15] GenomeInfoDb_1.26.7         IRanges_2.24.1             
[17] S4Vectors_0.28.1            BiocGenerics_0.36.1        
[19] MatrixGenerics_1.2.1        matrixStats_0.58.0         
[21] SeuratObject_4.0.0          Seurat_4.0.1               

loaded via a namespace (and not attached):
  [1] readxl_1.3.1              backports_1.2.1          
  [3] plyr_1.8.6                igraph_1.2.6             
  [5] lazyeval_0.2.2            splines_4.0.3            
  [7] BiocParallel_1.24.1       listenv_0.8.0            
  [9] scattermore_0.7           digest_0.6.27            
 [11] htmltools_0.5.1.1         viridis_0.6.0            
 [13] fansi_0.4.2               magrittr_2.0.1           
 [15] tensor_1.5                cluster_2.1.0            
 [17] ROCR_1.0-11               globals_0.14.0           
 [19] modelr_0.1.8              spatstat.sparse_2.0-0    
 [21] colorspace_2.0-0          rvest_1.0.0              
 [23] ggrepel_0.9.1             haven_2.4.1              
 [25] xfun_0.22                 crayon_1.4.1             
 [27] RCurl_1.98-1.3            jsonlite_1.7.2           
 [29] spatstat.data_2.1-0       survival_3.2-7           
 [31] zoo_1.8-9                 glue_1.4.2               
 [33] polyclip_1.10-0           gtable_0.3.0             
 [35] zlibbioc_1.36.0           XVector_0.30.0           
 [37] leiden_0.3.7              DelayedArray_0.16.3      
 [39] BiocSingular_1.6.0        future.apply_1.7.0       
 [41] abind_1.4-5               scales_1.1.1             
 [43] DBI_1.1.1                 miniUI_0.1.1.1           
 [45] Rcpp_1.0.6                viridisLite_0.4.0        
 [47] xtable_1.8-4              reticulate_1.19          
 [49] spatstat.core_2.1-2       rsvd_1.0.5               
 [51] htmlwidgets_1.5.3         httr_1.4.2               
 [53] RColorBrewer_1.1-2        ellipsis_0.3.2           
 [55] ica_1.0-2                 farver_2.1.0             
 [57] pkgconfig_2.0.3           scuttle_1.0.4            
 [59] dbplyr_2.1.1              sass_0.3.1               
 [61] uwot_0.1.10               deldir_0.2-10            
 [63] utf8_1.2.1                labeling_0.4.2           
 [65] tidyselect_1.1.1          rlang_0.4.10             
 [67] reshape2_1.4.4            later_1.2.0              
 [69] cellranger_1.1.0          munsell_0.5.0            
 [71] tools_4.0.3               cli_2.5.0                
 [73] generics_0.1.0            broom_0.7.6              
 [75] ggridges_0.5.3            evaluate_0.14            
 [77] fastmap_1.1.0             yaml_2.2.1               
 [79] goftest_1.2-2             fs_1.5.0                 
 [81] knitr_1.33                fitdistrplus_1.1-3       
 [83] RANN_2.6.1                pbapply_1.4-3            
 [85] future_1.21.0             nlme_3.1-149             
 [87] sparseMatrixStats_1.2.1   mime_0.10                
 [89] xml2_1.3.2                rstudioapi_0.13          
 [91] compiler_4.0.3            beeswarm_0.3.1           
 [93] plotly_4.9.3              png_0.1-7                
 [95] spatstat.utils_2.1-0      reprex_2.0.0             
 [97] bslib_0.2.4               stringi_1.5.3            
 [99] lattice_0.20-41           Matrix_1.2-18            
[101] vctrs_0.3.8               pillar_1.6.0             
[103] lifecycle_1.0.0           spatstat.geom_2.1-0      
[105] lmtest_0.9-38             jquerylib_0.1.4          
[107] RcppAnnoy_0.0.18          BiocNeighbors_1.8.2      
[109] data.table_1.14.0         cowplot_1.1.1            
[111] bitops_1.0-7              irlba_2.3.3              
[113] httpuv_1.6.0              patchwork_1.1.0.9000     
[115] R6_2.5.0                  promises_1.2.0.1         
[117] KernSmooth_2.23-17        gridExtra_2.3            
[119] vipor_0.4.5               parallelly_1.25.0        
[121] codetools_0.2-16          MASS_7.3-53              
[123] assertthat_0.2.1          withr_2.4.2              
[125] sctransform_0.3.2         GenomeInfoDbData_1.2.4   
[127] hms_1.0.0                 mgcv_1.8-33              
[129] grid_4.0.3                rpart_4.1-15             
[131] beachmat_2.6.4            rmarkdown_2.7            
[133] DelayedMatrixStats_1.12.3 Rtsne_0.15               
[135] lubridate_1.7.10          shiny_1.6.0              
[137] ggbeeswarm_0.6.0