Loading the R environment
library("DESeq")
library(reshape)
library(ggplot2)
Reading input
filename = "counttable1.tsv"
raw_counts = read.csv(filename, sep = "\t", header = T, row.names = 1, stringsAsFactors = F)
raw_counts_deseqct = newCountDataSet(raw_counts[, c(1, 2)], factor(c("can",
"can")))
Performing VST transformation
raw_counts_deseqct = newCountDataSet(raw_counts[, c(1, 2)], factor(c("can",
"can")))
raw_counts_deseqct = estimateSizeFactors(raw_counts_deseqct)
raw_counts_deseqctBlind = estimateDispersions(raw_counts_deseqct, method = "blind")
vsd = varianceStabilizingTransformation(raw_counts_deseqctBlind)
vst_raw_counts = exprs(vsd)
write.table(vst_raw_counts, file = "vst_transf_tabout.txt", sep = "\t", quote = F,
col.names = NA, row.names = TRUE)
head(vst_raw_counts)
## can_1 can_2
## FBgn0000003 10.095 10.155
## FBgn0000008 11.078 11.130
## FBgn0000014 9.834 9.921
## FBgn0000015 9.834 9.834
## FBgn0000017 12.450 12.406
## FBgn0000018 10.414 10.429
vst_raw_counts_melted = melt(vst_raw_counts)
ggplot(data = vst_raw_counts_melted, aes(value, colour = X2))
## Error: No layers in plot
+geom_density(alpha = 0.2)
## Error: invalid argument to unary operator
+labs(x = "VST of counts")
## Error: invalid argument to unary operator
for (sample in colnames(vst_raw_counts)) {
hist(vst_raw_counts[, colnames(vst_raw_counts) == sample], main = paste("Sample",
sample), breaks = 80, xlab = "VST counts")
}
Logging R environment
sessionInfo()
## R version 3.0.1 (2013-05-16)
## Platform: x86_64-pc-linux-gnu (64-bit)
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=C LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] parallel stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] ggplot2_0.9.3.1 reshape_0.8.4 plyr_1.8
## [4] DESeq_1.12.0 lattice_0.20-15 locfit_1.5-9.1
## [7] Biobase_2.20.1 BiocGenerics_0.6.0 markdown_0.6.1
## [10] knitr_1.2
##
## loaded via a namespace (and not attached):
## [1] annotate_1.38.0 AnnotationDbi_1.22.6 colorspace_1.2-2
## [4] DBI_0.2-7 dichromat_2.0-0 digest_0.6.3
## [7] evaluate_0.4.4 formatR_0.8 genefilter_1.42.0
## [10] geneplotter_1.38.0 grid_3.0.1 gtable_0.1.2
## [13] IRanges_1.18.2 labeling_0.2 MASS_7.3-27
## [16] munsell_0.4.2 proto_0.3-10 RColorBrewer_1.0-5
## [19] reshape2_1.2.2 RSQLite_0.11.4 scales_0.2.3
## [22] splines_3.0.1 stats4_3.0.1 stringr_0.6.2
## [25] survival_2.37-4 tools_3.0.1 XML_3.98-1.1
## [28] xtable_1.7-1