diff Dotplot_Release/Normalization_sigpreys.R @ 0:dfa3436beb67 draft

Uploaded
author bornea
date Fri, 29 Jan 2016 09:56:02 -0500
parents
children
line wrap: on
line diff
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/Dotplot_Release/Normalization_sigpreys.R	Fri Jan 29 09:56:02 2016 -0500
@@ -0,0 +1,46 @@
+#!/usr/bin/env Rscript
+
+#this programs normalizes a saint input file based on the spectral counts of "signficant" preys
+# that is, preys with an FDR <= the secondary cutoff as supplied to the dotplot script
+
+args <- commandArgs(trailingOnly = TRUE)
+
+d = read.delim(args[1], header=T, as.is=T)
+d <- d[d$BFDR <= as.numeric(args[2]),]
+
+baitn = 1
+curr_bait <- d$Bait[1]
+s <- vector()
+s[1] = 0
+for(i in 1:length(d$Bait)){
+	if(curr_bait != d$Bait[i]){
+		baitn <- baitn + 1
+		curr_bait <- d$Bait[i]
+		s[baitn] <- d$AvgSpec[i]
+	}
+	else{
+		s[baitn] <- s[baitn] + d$AvgSpec[i]
+	}
+}
+
+med.s = median(s)
+s = s / med.s
+
+d_n <- d
+baitn = 1
+curr_bait <- d_n$Bait[1]
+for(i in 1:length(d_n$Bait)){
+	if(curr_bait != d_n$Bait[i]){
+		baitn <- baitn + 1
+		curr_bait <- d_n$Bait[i]
+		d_n$AvgSpec[i] <- d_n$AvgSpec[i]/s[baitn]
+	}
+	else{
+		d_n$AvgSpec[i] <- d_n$AvgSpec[i]/s[baitn]
+	}
+}
+
+#print normalized data to file
+
+write.table(d_n, file = "norm_saint.txt", sep="\t", quote=F, row.names=F)
+