MBY2_CNV_stat.r 1.7 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445
  1. args<- commandArgs(TRUE)
  2. sam <- args[1]
  3. sam_data <- read.table(file=paste(sam,".cov",sep=""),header = TRUE,sep="\t",as.is=TRUE)
  4. sam_data <- sam_data[,c(1,5)]
  5. names(sam_data) <- c('contig',sam)
  6. sam_m <- mean(sam_data[,2])
  7. sam_sd<- sd(sam_data[,2])
  8. vals<-sam_data[,sam]
  9. vals=vals/sam_m
  10. vals[vals < 0.07] <- 0.07
  11. sam_data$Mean <- vals
  12. fivenum(sam_data[sam_data$contig=='Y',sam])
  13. y_mean <- mean(sam_data[sam_data$contig=='Y','Mean'])
  14. sex <- 'F'
  15. if(y_mean > 0.2){
  16. sex <- 'M'
  17. }
  18. pdf(file=paste(sam,sex,"QC_plot.pdf",sep="_"),width = 10)
  19. hist(sam_data[,sam],breaks=1000,xlim=c(0,1500),col='lightblue',
  20. xlab='Counts',main=paste(sam,' Counts/Bin. Mean= ',
  21. round(sam_m,0),' SD=',round(sam_sd,0),sep=''),freq = F)
  22. hist(sam_data[,'Mean'],breaks=1000,xlim=c(0,3),col='lightblue',
  23. xlab='Copy-Ratio',main=paste(sam,' Bin Copy-Ratio Distributions',sep=''),freq = F)
  24. hist(sam_data[sam_data$contig =='X' |sam_data$contig =='Y',sam],breaks=1000,xlim=c(0,1500),col='lightblue',
  25. xlab='Counts',main=paste(sam,' XY Counts/Bin. Mean= ',
  26. round(sam_m,0),' SD=',round(sam_sd,0),sep=''),freq = F)
  27. hist(sam_data[sam_data$contig =='X' |sam_data$contig =='Y','Mean'],breaks=1000,xlim=c(0,3),col='lightblue',
  28. xlab='Copy-Ratio',main=paste(sam,' XY Bin Copy-Ratio Distributions',sep=''),freq = F)
  29. boxplot(sam_data[,sam] ~ ordered(sam_data$contig,levels=c(1:22,'X','Y')),ylim=c(0,1500),
  30. main=paste(sam,'-Boxplot of Counts/Bin. Mean=',round(sam_m,0),
  31. 'SD/Mean=',round(sam_sd/sam_m,2),sep=' ') )
  32. lines(rep(sam_m,24),col='red')
  33. boxplot(sam_data[,'Mean'] ~ ordered(sam_data$contig,levels=c(1:22,'X','Y')),ylim=c(0,3),
  34. main=paste(sam,'-Boxplot of Copy-Ratio.',sep=' '))
  35. lines(rep(1,24),col='red')
  36. dev.off()