WES_CNV_stat.r 1.9 KB

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