figure_3_c_max.Rmd

``` {r setup, echo=FALSE, message=FALSE, include=FALSE, error=FALSE}
library(xtable)
library(ggplot2)
library(plyr)
library(BSgenome.Dmelanogaster.UCSC.dm3)
library(parallel)
library(rtracklayer)

# Output folder for this document
options(knitr.figure_dir = "figure_3_c_max_output")

source("shared_code/knitr_common.r")
source("shared_code/granges_common.r")
source("shared_code/samples.r")
source("shared_code/exo_metapeak.r")
source("shared_code/profiles_common.r")
```

# Figure 3C: Max ChIP-nexus

``` {r header_child, child="child_docs/header_child.Rmd"}
```

## Overview

We will compare the profiles of the following samples around the top 200 Max motifs:

``` {r max_samples_table, results="asis"}
max.df <- subset(samples.df, sample == "dmel_s2_max_chipnexus_01")
html_table(max.df)
```

## Motif

Max E-BOX motif: `CACGTG`

Window size centered at motif: 15 bp

``` {r calculate_profile_reads, include=FALSE}

reads.list <- cache("reads.list.rds", function() {
  max.motif <- "CACGTG"
  max.gr <- trim(filter_chrs(vmatchPattern(max.motif, Dmelanogaster, max.mismatch=0, fixed=FALSE)))
  max.gr <- max.gr[strand(max.gr) == "+"]

  checked_mclapply(max.df$sample, process_sample, max.gr, n=200, mc.cores=3)
})

```

## Average profile (top 200)

``` {r plots_per_sample, fig.cap="", fig.width=9, fig.height=6, dpi=100}
plots.list <- lapply(reads.list, build_plot)
nothing <- lapply(plots.list, print)
```

``` {r zoomed_plot, warning=FALSE, fig.cap="", fig.width=9, fig.height=6, dpi=100}
reads.df <- reads.list[[1]]

motif.box <- data.frame(xmin=0, 
                        xmax=reads.df$motif_width[1]-1,
                        ymin=-Inf,
                        ymax=Inf)
                          
g <- ggplot(reads.df, aes(x=tss_distance, y=reads, color=strand)) +
     geom_line(size=1.2) +
     geom_vline(xintercept=0:5, color="gray50") +
     scale_colour_manual(name="Strand", values=c("+"="red", "-"="darkblue")) +
     geom_rect(show_guide=FALSE, inherit.aes=FALSE, data=motif.box, 
               aes(xmin=xmin, xmax=xmax, ymin=ymin, ymax=ymax), 
               alpha=0.25, fill="gray80") +
     theme_bw() +
     scale_x_continuous(breaks=c(-15, -10, -8, -5, 0:5, 10, 13, 15, 20), limits=c(-15, 21)) +
     labs(x="Distance to motif left edge", y="Average ChIP-nexus reads", title="CACGTG")
g
```

``` {r save_zoomed_pdf, include=FALSE}
pdf(figure_path("max_CACGTG_profile.pdf"), width=9, height=6)
print(g)
dev.off()
```

## Heatmap

``` {r collect_heatmap_reads, include=FALSE}

beds <- list.files(figure_path(), "bed")
names(beds) <- gsub("_top_200_motifs.bed", "", beds)

sample.list <- mclapply(beds, heatmap_reads, mc.cores=4)
```

``` {r plot_heatmaps, fig.cap="", fig.width=6, fig.heigth=6, dpi=100}
nothing <- lapply(names(sample.list), function(n) { draw_exo_heatmap(sample.list[[n]], n, 75:125)})
```

``` {r create_heatmap_pdf, include=FALSE}
pdf(figure_path("max_heatmaps.pdf"), width=6, height=6)
nothing <- lapply(names(sample.list), function(n) { draw_exo_heatmap(sample.list[[n]], n, 75:125)})
dev.off()
```


``` {r session_info_child, child="child_docs/session_info_child.Rmd"}
```