Commit 838fe7f8 authored by Holger Brandl's avatar Holger Brandl

fixed ngs-plot

parent 2109980a
......@@ -155,7 +155,7 @@ degs <- deResults %>% filter(is_hit)
ggplot(degs, aes(paste(sample_1, "vs", sample_2))) + geom_bar() +xlab(NULL) + ylab("# DBGs") +ggtitle("DBG count summary") + coord_flip()
ggplot(degs, aes(paste(sample_1, "vs", sample_2), fill=log2FoldChange>1)) + geom_bar(position="dodge") +xlab(NULL) + ylab("# DBGs") +ggtitle("DBG count summary by direction of expression") + coord_flip()
ggplot(degs, aes(paste(sample_1, "vs", sample_2), fill=sample_1_overex)) + geom_bar(position="dodge") +xlab(NULL) + ylab("# DBGs") +ggtitle("DBG count summary by direction of expression") + coord_flip()
## just needed to restore environment easily
save(degs, file=".degs.RData")
......@@ -165,6 +165,12 @@ save(degs, file=".degs.RData")
#res %>% as.df() %>% ggplot(aes(padj)) + geom_histogram()
## todo add gene info
degs %>% inner_join(geneInfo) %>%write.delim(file="degs.txt")
# degs <- read.delim("degs.txt")
#' [degs](degs.txt)
plotMA(res, main="DESeq2", ylim=c(-2,2))
......
......@@ -40,6 +40,10 @@ splitProtTime <- . %>% separate(sample, c("protein", "timpoint"), sep="_", remov
peaks %<>% mutate(peak_width=end_position-start_position)
save(peaks, file="peaks.RData")
# peaks <- local(get(load("peaks.RData")))
peaks %>%
# sample_frac(0.3) %>%
ggplot(aes(peak_width)) +
......@@ -91,6 +95,8 @@ peaks %>% sample_frac(0.1) %>%
arrange(-score) %>%
mutate(x=row_number()/n()) %>%
# note: max: score ~ 10*qvalue
# sample_frac(0.3) %>%
ggplot(aes(x, qvalue, color=sample)) +
geom_line() +
......@@ -98,16 +104,61 @@ peaks %>% sample_frac(0.1) %>%
xlim(0, 0.2) + ylim(0,50)
#
## chippeak anno
########################################################################################################################
#' # ChippeakAnno
# see [docs](http://www.bioconductor.org/packages/release/bioc/vignettes/ChIPpeakAnno/inst/doc/ChIPpeakAnno.pdf)
#+ eval=FALSE, echo=FALSE
## DEBUG start
if(F){
#http://www.bioconductor.org/packages/release/bioc/vignettes/ChIPpeakAnno/inst/doc/ChIPpeakAnno.pdf
require(ChIPpeakAnno)
data(TSS.zebrafish.Zv9)
## https://www.biostars.org/p/115101/
mydata.ranged<-BED2RangedData("hmedip_allpeaks.bed",header=FALSE)
## http://master.bioconductor.org/help/course-materials/2009/SeattleNov09/IRanges/IRangesOverview.R
#http://www.bioconductor.org/packages/release/bioc/vignettes/ChIPpeakAnno/inst/doc/ChIPpeakAnno.pdf
peakSet <- peaks %>% sample_n(10000) %>% filter(sample=="H3HA_Dome", peak_type=="narrow_peaks")
peakRanges <-with(peakSet, RangedData(ranges = IRanges(start_position, end_position), strand = ".", space = chromosome_name))
### other feature types
#Anno_exon = getAnnotation(mart, featureType=c("Exon"))
#Anno_TSS = getAnnotation(mart, featureType=c("TSS"))
#Anno_miRNA = getAnnotation(mart, featureType=c("miRNA"))
#annotatedPeak_exon = annotatePeakInBatch(pr, AnnotationData=Anno_exon)
annotatedPeak = annotatePeakInBatch (peakRanges, AnnotationData = TSS.zebrafish.Zv9)
annotatedPeak %>% as.df()%>% ggplot(aes(insideFeature)) + geom_bar()
}
## DEBUG end
#+
require(ChIPpeakAnno)
data(TSS.zebrafish.Zv9)
## apply csa to all sets
annotatePeaks <- function(somePeaks, subsample=min(nrow(peaks), 10000)){
peakRanges <- somePeaks %>%
sample_n(subsample) %$%
RangedData(ranges = IRanges(start_position, end_position), strand = ".", space = chromosome_name)
annotatePeakInBatch (peakRanges, AnnotationData = TSS.zebrafish.Zv9)
annotatedPeak %>% as.df()
}
## DEBUG end
\ No newline at end of file
annoPeaks <- peaks %>% group_by(sample, peak_type) %>% do(annotatePeaks(.))
save(annoPeaks, file="annoPeaks.RData")
# annoPeaks <- local(get(load("annoPeaks.RData")))
with(peaks, as.data.frame(table(sample, peak_type)))
peaks %$% as.data.frame(table(sample, peak_type)) %>% head
annoPeaks %>% head(20) %>% kable()
annoPeaks %>% ungroup %$%
as.data.frame(table(sample, peak_type, insideFeature)) %>% head
ggplot(aes(sample, fill=insideFeature)) + geom_bar() + facet_wrap(~peak_type)
\ No newline at end of file
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment