diff --git a/README.Rmd b/README.Rmd
index a47a037..08de0dd 100644
--- a/README.Rmd
+++ b/README.Rmd
@@ -57,11 +57,8 @@ First we'll also need to load a couple dependencies & resolve a function conflic
 
 ```{r libraries, results='hide', message=FALSE}
 library(dplyr)
-library(scater)
 library(scLANE)
 library(ggplot2)
-select <- dplyr::select
-filter <- dplyr::filter
 ```
 
 ## Input data 
@@ -75,14 +72,23 @@ sim_data <- readRDS(url("https://zenodo.org/records/8433077/files/scLANE_sim_dat
 The PCA embeddings show us a pretty simple trajectory that's strongly correlated with the first principal component. 
 
 ```{r plot-sims-pt, results='hold'}
-plotPCA(sim_data, colour_by = "cell_time_normed") + 
+data.frame(sim_data@int_colData$reducedDims@listData$PCA[, 1:2]) %>% 
+  mutate(pseudotime = sim_data$cell_time_normed) %>% 
+  ggplot(aes(x = PC1, y = PC2, color = pseudotime)) + 
+  geom_point(size = 2, alpha = 0.75, stroke = 0) + 
+  scale_color_gradientn(colors = viridisLite::plasma(n = 20)) + 
+  labs(x = "PC 1", y = "PC 2", color = "Pseudotime") + 
   theme_scLANE(umap = TRUE)
 ```
 
 We also see that the data are not clustered by subject, which indicates that gene dynamics are mostly homogeneous across subjects. 
 
 ```{r plot-sims-subj, results='hold'}
-plotPCA(sim_data, colour_by = "subject") + 
+data.frame(sim_data@int_colData$reducedDims@listData$PCA[, 1:2]) %>% 
+  mutate(subject = sim_data$subject) %>% 
+  ggplot(aes(x = PC1, y = PC2, color = subject)) + 
+  geom_point(size = 2, alpha = 0.75, stroke = 0) + 
+  labs(x = "PC 1", y = "PC 2", color = "Subject ID") + 
   theme_scLANE(umap = TRUE)
 ```