Skip to content

Latest commit

 

History

History
54 lines (39 loc) · 1.86 KB

002.Figure 1A (ATAC_PCA).md

File metadata and controls

54 lines (39 loc) · 1.86 KB

PCA

## calculate mapped reads in each region of 11 chromosomes for opium poppy
bedtools makewindows -g poppy.genome -w 20000 > poppy.windows
multiBamSummary BED-file  --BED ./chr/poppy.windows --bamfiles ./*.bam  --numberOfProcessors 20 -out ./chr/readCounts.npz --outRawCounts ./chr/readCounts.tab

library(FactoMineR)
library(corrplot)
library(factoextra)
atac_seq <- read.csv("D:/001poppy_atac_new_genome/ATAC_peak/atac_peak_cor/readCounts.tab",sep = "\t",stringsAsFactors = F)


atac_seq <- atac_seq[,-1:-3]


colnames(atac_seq) <- c("Capsule_1","Capsule_2","Capsule_3",
                        "Stem_1","Stem_2","Stem_3",
                        "Fine Root_1","Fine Root_2","Fine Root_3",
                        "Leaf_1","Leaf_2","Leaf_3",
                        "Tap Root_1","Tap Root_2","Tap Root_3",
                        "Petal_1","Petal_2","Petal_3")

plot_atac_seq <- as.data.frame(t(atac_seq))
dim(plot_atac_seq)
plot_atac_seq[,108632] <- c("Capsule","Capsule","Capsule",
                            "Stem","Stem","Stem",
                            "Fine Root","Fine Root","Fine Root",
                            "Leaf","Leaf","Leaf",
                            "Tap Root","Tap Root","Tap Root",
                            "Petal","Petal","Petal")


colnames(plot_atac_seq)[108632] <- "Tissue"


atac_seq_pca <- FactoMineR::PCA(t(atac_seq),scale.unit = FALSE, graph = FALSE,ncp=10)

fviz_eig(atac_seq_pca,addlabels = T)
fviz_pca_ind(atac_seq_pca,axes = c(1,2),labelsize =6, 
             habillage = factor(plot_atac_seq$Tissue),
             palette = c("#E41A1C", "#377EB8" ,"#4DAF4A" ,"#984EA3" ,"#FF7F00"  ,"#A65628"),
             invisible='quali'
) +
  ggforce::geom_mark_ellipse(aes(fill = Groups,
                                 color = Groups)) +
  theme(legend.position = 'bottom') +
  coord_equal()