comparison RScript.r @ 3:d27d745d0556 draft

Uploaded
author davidvanzessen
date Mon, 14 Oct 2013 07:20:38 -0400
parents aac19ac7cb21
children 10cfa5e9186e
comparison
equal deleted inserted replaced
2:aac19ac7cb21 3:d27d745d0556
109 109
110 plotVD <- function(dat){ 110 plotVD <- function(dat){
111 ggplot() + 111 ggplot() +
112 geom_tile(data=dat, aes(x=factor(Top.D.Gene), y=factor(Top.V.Gene), fill=log)) + 112 geom_tile(data=dat, aes(x=factor(Top.D.Gene), y=factor(Top.V.Gene), fill=log)) +
113 theme(axis.text.x = element_text(angle = 90, hjust = 1)) + 113 theme(axis.text.x = element_text(angle = 90, hjust = 1)) +
114 scale_fill_gradient(low="gold", high="blue") + 114 scale_fill_gradient(low="gold", high="blue", na.value="white") +
115 ggtitle(unique(dat$Sample)) + 115 ggtitle(unique(dat$Sample)) +
116 xlab("D genes") + 116 xlab("D genes") +
117 ylab("V Genes") 117 ylab("V Genes")
118 } 118 }
119 119
122 cartegianProductVD = expand.grid(Top.V.Gene = Vchain$v.name, Top.D.Gene = Dchain$v.name, Sample = unique(test$Sample)) 122 cartegianProductVD = expand.grid(Top.V.Gene = Vchain$v.name, Top.D.Gene = Dchain$v.name, Sample = unique(test$Sample))
123 123
124 completeVD = merge(VandDCount, cartegianProductVD, all.y=TRUE) 124 completeVD = merge(VandDCount, cartegianProductVD, all.y=TRUE)
125 completeVD$Length = as.numeric(completeVD$Length) 125 completeVD$Length = as.numeric(completeVD$Length)
126 completeVD$log = log(completeVD$Length) 126 completeVD$log = log(completeVD$Length)
127 completeVD$log[is.na(completeVD$log)] = 0 127 #completeVD$log[is.na(completeVD$log)] = 0
128 l = split(completeVD, f=completeVD[,"Sample"]) 128 l = split(completeVD, f=completeVD[,"Sample"])
129 png("HeatmapVD%d.png", width=150+(15*length(Dchain$v.name)), height=100+(15*length(Vchain$v.name))) 129 png("HeatmapVD%d.png", width=150+(15*length(Dchain$v.name)), height=100+(15*length(Vchain$v.name)))
130 lapply(l, FUN=plotVD) 130 lapply(l, FUN=plotVD)
131 dev.off() 131 dev.off()
132 132
133 133
134 plotVJ <- function(dat){ 134 plotVJ <- function(dat){
135 ggplot() + 135 ggplot() +
136 geom_tile(data=dat, aes(x=factor(Top.J.Gene), y=factor(Top.V.Gene), fill=log)) + 136 geom_tile(data=dat, aes(x=factor(Top.J.Gene), y=factor(Top.V.Gene), fill=log)) +
137 theme(axis.text.x = element_text(angle = 90, hjust = 1)) + 137 theme(axis.text.x = element_text(angle = 90, hjust = 1)) +
138 scale_fill_gradient(low="gold", high="blue") + 138 scale_fill_gradient(low="gold", high="blue", na.value="white") +
139 ggtitle(unique(dat$Sample)) + 139 ggtitle(unique(dat$Sample)) +
140 xlab("J genes") + 140 xlab("J genes") +
141 ylab("V Genes") 141 ylab("V Genes")
142 } 142 }
143 143
145 cartegianProductVJ = expand.grid(Top.V.Gene = Vchain$v.name, Top.J.Gene = Jchain$v.name, Sample = unique(test$Sample)) 145 cartegianProductVJ = expand.grid(Top.V.Gene = Vchain$v.name, Top.J.Gene = Jchain$v.name, Sample = unique(test$Sample))
146 146
147 completeVJ = merge(VandJCount, cartegianProductVJ, all.y=TRUE) 147 completeVJ = merge(VandJCount, cartegianProductVJ, all.y=TRUE)
148 completeVJ$Length = as.numeric(completeVJ$Length) 148 completeVJ$Length = as.numeric(completeVJ$Length)
149 completeVJ$log = log(completeVJ$Length) 149 completeVJ$log = log(completeVJ$Length)
150 completeVJ$log[is.na(completeVJ$log)] = 0 150 #completeVJ$log[is.na(completeVJ$log)] = 0
151 l = split(completeVJ, f=completeVJ[,"Sample"]) 151 l = split(completeVJ, f=completeVJ[,"Sample"])
152 png("HeatmapVJ%d.png", width=150+(15*length(Jchain$v.name)), height=100+(15*length(Vchain$v.name))) 152 png("HeatmapVJ%d.png", width=150+(15*length(Jchain$v.name)), height=100+(15*length(Vchain$v.name)))
153 lapply(l, FUN=plotVJ) 153 lapply(l, FUN=plotVJ)
154 dev.off() 154 dev.off()
155 155
156 plotDJ <- function(dat){ 156 plotDJ <- function(dat){
157 ggplot() + 157 ggplot() +
158 geom_tile(data=dat, aes(x=factor(Top.J.Gene), y=factor(Top.D.Gene), fill=log)) + 158 geom_tile(data=dat, aes(x=factor(Top.J.Gene), y=factor(Top.D.Gene), fill=log)) +
159 theme(axis.text.x = element_text(angle = 90, hjust = 1)) + 159 theme(axis.text.x = element_text(angle = 90, hjust = 1)) +
160 scale_fill_gradient(low="gold", high="blue") + 160 scale_fill_gradient(low="gold", high="blue", na.value="white") +
161 ggtitle(unique(dat$Sample)) + 161 ggtitle(unique(dat$Sample)) +
162 xlab("J genes") + 162 xlab("J genes") +
163 ylab("D Genes") 163 ylab("D Genes")
164 } 164 }
165 165
167 cartegianProductDJ = expand.grid(Top.D.Gene = Dchain$v.name, Top.J.Gene = Jchain$v.name, Sample = unique(test$Sample)) 167 cartegianProductDJ = expand.grid(Top.D.Gene = Dchain$v.name, Top.J.Gene = Jchain$v.name, Sample = unique(test$Sample))
168 168
169 completeDJ = merge(DandJCount, cartegianProductDJ, all.y=TRUE) 169 completeDJ = merge(DandJCount, cartegianProductDJ, all.y=TRUE)
170 completeDJ$Length = as.numeric(completeDJ$Length) 170 completeDJ$Length = as.numeric(completeDJ$Length)
171 completeDJ$log = log(completeDJ$Length) 171 completeDJ$log = log(completeDJ$Length)
172 completeDJ$log[is.na(completeDJ$log)] = 0 172 #completeDJ$log[is.na(completeDJ$log)] = 0
173 l = split(completeDJ, f=completeDJ[,"Sample"]) 173 l = split(completeDJ, f=completeDJ[,"Sample"])
174 png("HeatmapDJ%d.png", width=150+(15*length(Jchain$v.name)), height=100+(15*length(Dchain$v.name))) 174 png("HeatmapDJ%d.png", width=150+(15*length(Jchain$v.name)), height=100+(15*length(Dchain$v.name)))
175 lapply(l, FUN=plotDJ) 175 lapply(l, FUN=plotDJ)
176 dev.off() 176 dev.off()
177 177