DEseq2_PCA_manually.md

## DEseq2 built-in function
plotPCA(vsd.fast, intgroup=c("subtype"))

##SVD to get PCs mannually
X<- assay(vsd.fast)

## center X
X<- t(scale(t(X),center=TRUE,scale=FALSE))

sv<- svd(t(X))
U<- sv$u
V<- sv$v
D<- sv$d
Z<- t(X)%*%V

variance_explained_each_PC<- function(x){
        var.list= list()
        varex = 0
        cumvar = 0
        denom = sum(x^2)
        for(i in 1:length(x)){
                varex[i] = x[i]^2/denom
                cumvar[i] = sum(x[1:i]^2)/denom
        }
        var.list$varex<- varex
        var.list$cumvar<- cumvar
        var.list
}

### screen plot 
screen.plot<- function(var.list){
        par(mfrow=c(1,2))
        plot(1:length(var.list$varex), var.list$varex *100,type="h",
             lwd=2,xlab="PC",ylab="% Variance Explained")
        plot(1:length(var.list$cumvar),var.list$cumvar,type="h",
             lwd=2,xlab="PC",ylab="Cummulative Variance Explained")
        
}

screen.plot(variance_explained_each_PC(D))

pc_dat<- data.frame(cell.line = coldata$subtype, PC1 = Z[,1], PC2 = Z[,2])

## make figure with ggplot2
library(ggplot2)
library(ggthemes)
ggplot(pc_dat, aes(x=PC1, y=PC2, col=cell.line, shape = cell.line)) + 
        geom_point(size = 3) + 
        scale_shape_manual(values=c(19, 4, 8, 3)) + 
        theme_bw(base_size = 14) +
        theme(panel.border = element_blank(), 
              panel.grid.major = element_blank(),
              panel.grid.minor = element_blank(), 
              axis.line.x = element_line(color="black", size = 0.6),
              axis.line.y = element_line(color="black", size = 0.6))