added expression filter clustering

dge_workflow/dge_analysis.R

@@ -70,9 +70,9 @@ runInfo(cuff)
 replicateInfo <- replicates(cuff) %>% mutate(file=basename(file)) %>% select(-c(total_mass, norm_mass, internal_scale, external_scale))
 replicateInfo
 
-write.delim(replicateInfo, file="replicateInfo.txt")
-# replicateInfo <- read.delim("replicateInfo.txt")
-#'  [replicateInfo](replicateInfo.txt)
+write.delim(replicateInfo, file="replicate_info.txt")
+# replicateInfo <- read.delim("replicate_info.txt")
+#'  [replicateInfo](replicate_info.txt)
 
 #+ eval=skipQC
 #hello hidden field
@@ -117,27 +117,20 @@ repFpkmMatrix(genes(cuff)) %>%
 # gene_fpkms <- read.delim("gene_fpkms_by_replicate.txt")
 #' [Annotated Expresssion Table](gene_fpkms_by_replicate.txt)
 
-fpkm(genes(cuff)) %>%
-    dplyr::rename(ensembl_gene_id=gene_id) %>%
-    merge(geneInfo, by="ensembl_gene_id", all.x=T) %>%
-    write.delim(file="geneFpkm.txt")
-
-
 ## also export count tables
-
-count(genes(cuff)) %>%
-    dplyr::rename(ensembl_gene_id=gene_id) %>%
+countMatrix(genes(cuff)) %>%
+    rownames2column("ensembl_gene_id") %>%
     merge(geneInfo, by="ensembl_gene_id", all.x=T) %>%
     write.delim(file="gene_counts.txt")
-# gene_counts <- read.delim("gene_counts.txt")
-#'  [Gene Fragment Counts](gene_counts.txt)
 
-countMatrix(genes(cuff)) %>%
+#'  [Gene Counts](gene_counts.txt)
+
+repCountMatrix(genes(cuff)) %>%
     rownames2column("ensembl_gene_id") %>%
     merge(geneInfo, by="ensembl_gene_id", all.x=T) %>%
-    write.delim(file="replicateCounts.txt")
+    write.delim(file="gene_counts_by_replicate.txt")
 # repCounts <- read.delim("replicateCounts.txt")
-#'  [Raw Counts by Replicate](replicateCounts.txt)
+#'  [Gene Counts by Replicate](replicateCounts.txt)
 
 ##ggplot(geneFpkm, aes(fpkm)) + geom_histogram() + scale_x_log10()
 
@@ -205,6 +198,25 @@ csDendro(genes(cuff),replicates=T)
 #dev.off()
 
 
+#' Also cluster expression filtered data
+exprGenes <- repFpkmMatrix(genes(cuff)) %>%
+    rownames2column("ensembl_gene_id") %>%
+    melt(value.name="fpkm") %>%
+    group_by(ensembl_gene_id) %>%
+    ## apply log trafo
+    mutate(log_fpkm=log10(fpkm+1)) %>%
+    ## apply expression filter
+    filter(max(log_fpkm)>1) %>%
+    dcast( ensembl_gene_id ~ variable, value.var="log_fpkm") %>%
+    column2rownames("ensembl_gene_id")
+
+res<-JSdist(makeprobs(exprGenes)) %>% hclust() %>% as.dendrogram()
+
+par(mar=c(10,4,4,4))
+plot(res, main="Replicate Clustering of Expressed Genes")
+
+
+
 
 #######################################################################################################################
 #' ## Extract lists of differentially expressed genes