Mercurial > repos > davidvanzessen > mutation_analysis

comparison aa_histogram.r @ 56:a4317b006d70 draft

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Uploaded
author davidvanzessen
date Tue, 01 Mar 2016 05:06:47 -0500
parents d4e72eeea640
children cb66d6dd1e66
comparison
equal deleted inserted replaced
55:0d5add1a9800 56:a4317b006d70
3 args <- commandArgs(trailingOnly = TRUE) 3 args <- commandArgs(trailingOnly = TRUE)
4 4
5 input = args[1] 5 input = args[1]
6 outfile = args[2] 6 outfile = args[2]
7 7
8 print("---------------- read input ----------------")
9
8 dat = read.table(input, sep="\t", fill=T, header=T, quote="") 10 dat = read.table(input, sep="\t", fill=T, header=T, quote="")
9 11
10 12 print("---------------- as numeric ----------------")
11 13
12 mutations.at.position = as.numeric(dat[1,]) 14 mutations.at.position = as.numeric(dat[1,])
13 aa.at.position = as.numeric(dat[2,]) 15 aa.at.position = as.numeric(dat[2,])
14 16
17 print("---------------- freq data.frame ----------------")
18
15 dat_freq = mutations.at.position / aa.at.position 19 dat_freq = mutations.at.position / aa.at.position
16 dat_dt = data.frame(i=1:length(dat_freq), freq=dat_freq) 20 dat_dt = data.frame(i=1:length(dat_freq), freq=dat_freq)
21
22 print("---------------- weird stuff ----------------")
17 23
18 options(width=220) 24 options(width=220)
19 25
20 print(dat[,20:40]) 26 print(dat[,20:40])
21 27
22 print(dat_dt) #need this or it will fail???? 28 print(dat_dt) #need this or it will fail????
29
30 print("---------------- plot ----------------")
23 31
24 m = ggplot(dat_dt, aes(x=i, y=freq)) + theme(axis.text.x = element_text(angle = 90, hjust = 1)) 32 m = ggplot(dat_dt, aes(x=i, y=freq)) + theme(axis.text.x = element_text(angle = 90, hjust = 1))
25 m = m + geom_histogram(stat="identity", colour = "black", fill = "darkgrey", alpha=0.8) + scale_x_continuous(breaks=1:length(dat_freq), labels=1:length(dat_freq)) 33 m = m + geom_histogram(stat="identity", colour = "black", fill = "darkgrey", alpha=0.8) + scale_x_continuous(breaks=1:length(dat_freq), labels=1:length(dat_freq))
26 m = m + annotate("segment", x = 0.5, y = -0.05, xend=26.5, yend=-0.05, colour="darkgreen", size=1) + annotate("text", x = 13, y = -0.1, label="FR1") 34 m = m + annotate("segment", x = 0.5, y = -0.05, xend=26.5, yend=-0.05, colour="darkgreen", size=1) + annotate("text", x = 13, y = -0.1, label="FR1")
27 m = m + annotate("segment", x = 26.5, y = -0.07, xend=38.5, yend=-0.07, colour="darkblue", size=1) + annotate("text", x = 32.5, y = -0.15, label="CDR1") 35 m = m + annotate("segment", x = 26.5, y = -0.07, xend=38.5, yend=-0.07, colour="darkblue", size=1) + annotate("text", x = 32.5, y = -0.15, label="CDR1")