Visulz_bulkATAC_PCA.R

# MESSAGE -----------------------------------------------------------------
#
# author: Yulin Lyu
# email: lvyulin@pku.edu.cn
#
# require: R whatever
#
# ---

# * 1. Load packages ------------------------------------------------------

setwd("exampleData/ATAC")

# grammar
library(tidyverse)
library(magrittr)
library(glue)
library(data.table)

# analysis
library(irlba)

# graphics
library(ggplot2)
library(ggrepel)
library(ggsci)
library(scales)

# * 2. Load data ----------------------------------------------------------

peakMtx <- readRDS("peakMtx.rds")

# * 3. Plot ---------------------------------------------------------------

dir.create("graphics")

PCAdata <- prcomp_irlba(peakMtx, n = 3, scale. = T)
PCprop <- (PCAdata$sdev)^2 %>% {round(. / sum(.), 3) * 100}

plotData <- as.data.table(PCAdata$rotation)
plotData[, id := colnames(peakMtx)][, type := rep(c("F", "P", "XF"), c(2, 3, 3))][]

ggplot(plotData, aes(x = PC1, y = PC2)) +
  geom_point(aes(color = type), show.legend = F) +
  geom_text_repel(aes(label = id)) +
  scale_color_d3() +
  labs(
    x = str_c("PC1 (", PCprop[1], "%)"),
    y = str_c("PC2 (", PCprop[2], "%)")) +
  theme(
    aspect.ratio = 1,
    panel.background = element_blank(),
    panel.grid = element_blank(),
    panel.border = element_rect(fill = NA)
  )

ggsave("graphics/PCA.png", width = 4, height = 4)