Contents
Set-up
library(BiocStyle)
library(knitr)
options(digits=3)
options(width=90)
url <- "http://www.ncbi.nlm.nih.gov/geo/download/?acc=GSE63310&format=file"
utils::download.file(url, destfile="GSE63310_RAW.tar", mode="wb")
utils::untar("GSE63310_RAW.tar", exdir = ".")
files <- c("GSM1545535_10_6_5_11.txt", "GSM1545536_9_6_5_11.txt", "GSM1545538_purep53.txt",
"GSM1545539_JMS8-2.txt", "GSM1545540_JMS8-3.txt", "GSM1545541_JMS8-4.txt",
"GSM1545542_JMS8-5.txt", "GSM1545544_JMS9-P7c.txt", "GSM1545545_JMS9-P8c.txt")
for(i in paste(files, ".gz", sep=""))
R.utils::gunzip(i, overwrite=TRUE)
Note that this workflow uses the developmental version of the Glimma package, so be ensure to install as follows before running.
source("http://www.bioconductor.org/biocLite.R")
## Bioconductor version 3.4 (BiocInstaller 1.23.4), ?biocLite for help
biocLite("Glimma", suppressUpdates=TRUE)
## BioC_mirror: https://bioconductor.org
## Using Bioconductor 3.4 (BiocInstaller 1.23.4), R 3.3.0 (2016-05-03).
## Installing package(s) 'Glimma'
##
## The downloaded binary packages are in
## /var/folders/lv/7wb15fhx4cx2s3sxn2wxtzh000010v/T//RtmpqoUyB0/downloaded_packages
Data packaging
files <- c("GSM1545535_10_6_5_11.txt", "GSM1545536_9_6_5_11.txt",
"GSM1545538_purep53.txt", "GSM1545539_JMS8-2.txt",
"GSM1545540_JMS8-3.txt", "GSM1545541_JMS8-4.txt",
"GSM1545542_JMS8-5.txt", "GSM1545544_JMS9-P7c.txt",
"GSM1545545_JMS9-P8c.txt")
read.delim(files[1], nrow=5)
## EntrezID GeneLength Count
## 1 497097 3634 1
## 2 100503874 3259 0
## 3 100038431 1634 0
## 4 19888 9747 0
## 5 20671 3130 1
library(limma)
library(edgeR)
x <- readDGE(files, columns=c(1,3))
class(x)
## [1] "DGEList"
## attr(,"package")
## [1] "edgeR"
dim(x)
## [1] 27179 9
samplenames <- substring(colnames(x), 12, nchar(colnames(x)))
samplenames
## [1] "10_6_5_11" "9_6_5_11" "purep53" "JMS8-2" "JMS8-3" "JMS8-4" "JMS8-5"
## [8] "JMS9-P7c" "JMS9-P8c"
colnames(x) <- samplenames
group <- as.factor(c("LP", "ML", "Basal", "Basal", "ML", "LP",
"Basal", "ML", "LP"))
x$samples$group <- group
lane <- as.factor(rep(c("L004","L006","L008"), c(3,4,2)))
x$samples$lane <- lane
x$samples
## files group lib.size norm.factors lane
## 10_6_5_11 GSM1545535_10_6_5_11.txt LP 32863052 1 L004
## 9_6_5_11 GSM1545536_9_6_5_11.txt ML 35335491 1 L004
## purep53 GSM1545538_purep53.txt Basal 57160817 1 L004
## JMS8-2 GSM1545539_JMS8-2.txt Basal 51368625 1 L006
## JMS8-3 GSM1545540_JMS8-3.txt ML 75795034 1 L006
## JMS8-4 GSM1545541_JMS8-4.txt LP 60517657 1 L006
## JMS8-5 GSM1545542_JMS8-5.txt Basal 55086324 1 L006
## JMS9-P7c GSM1545544_JMS9-P7c.txt ML 21311068 1 L008
## JMS9-P8c GSM1545545_JMS9-P8c.txt LP 19958838 1 L008
library(Mus.musculus)
geneid <- rownames(x)
genes <- select(Mus.musculus, keys=geneid, columns=c("SYMBOL", "TXCHROM"),
keytype="ENTREZID")
head(genes)
## ENTREZID SYMBOL TXCHROM
## 1 497097 Xkr4 chr1
## 2 100503874 Gm19938 <NA>
## 3 100038431 Gm10568 <NA>
## 4 19888 Rp1 chr1
## 5 20671 Sox17 chr1
## 6 27395 Mrpl15 chr1
dup <- genes$ENTREZID[duplicated(genes$ENTREZID)]
genes[genes$ENTREZID %in% dup,][1:10,]
## ENTREZID SYMBOL TXCHROM
## 5360 100316809 Mir1906-1 chr12
## 5361 100316809 Mir1906-1 chrX
## 7128 218963 Gm1821 chr11
## 7129 218963 Gm1821 chr14
## 9564 12228 Btg3 chr16
## 9565 12228 Btg3 chr17
## 11351 433182 Eno1b chr4
## 11352 433182 Eno1b chr18
## 11544 100217457 Snord58b chr14
## 11545 100217457 Snord58b chr18
mat <- match(geneid, genes$ENTREZID)
genes <- genes[mat,]
genes[genes$ENTREZID %in% dup,][1:5,]
## ENTREZID SYMBOL TXCHROM
## 5360 100316809 Mir1906-1 chr12
## 7128 218963 Gm1821 chr11
## 9564 12228 Btg3 chr16
## 11351 433182 Eno1b chr4
## 11544 100217457 Snord58b chr14
x$genes <- genes
x
## An object of class "DGEList"
## $samples
## files group lib.size norm.factors lane
## 10_6_5_11 GSM1545535_10_6_5_11.txt LP 32863052 1 L004
## 9_6_5_11 GSM1545536_9_6_5_11.txt ML 35335491 1 L004
## purep53 GSM1545538_purep53.txt Basal 57160817 1 L004
## JMS8-2 GSM1545539_JMS8-2.txt Basal 51368625 1 L006
## JMS8-3 GSM1545540_JMS8-3.txt ML 75795034 1 L006
## JMS8-4 GSM1545541_JMS8-4.txt LP 60517657 1 L006
## JMS8-5 GSM1545542_JMS8-5.txt Basal 55086324 1 L006
## JMS9-P7c GSM1545544_JMS9-P7c.txt ML 21311068 1 L008
## JMS9-P8c GSM1545545_JMS9-P8c.txt LP 19958838 1 L008
##
## $counts
## Samples
## Tags 10_6_5_11 9_6_5_11 purep53 JMS8-2 JMS8-3 JMS8-4 JMS8-5 JMS9-P7c JMS9-P8c
## 497097 1 2 342 526 3 3 535 2 0
## 100503874 0 0 5 6 0 0 5 0 0
## 100038431 0 0 0 0 0 0 1 0 0
## 19888 0 1 0 0 17 2 0 1 0
## 20671 1 1 76 40 33 14 98 18 8
## 27174 more rows ...
##
## $genes
## ENTREZID SYMBOL TXCHROM
## 1 497097 Xkr4 chr1
## 2 100503874 Gm19938 <NA>
## 3 100038431 Gm10568 <NA>
## 4 19888 Rp1 chr1
## 5 20671 Sox17 chr1
## 27174 more rows ...
Data pre-processing
cpm <- cpm(x)
lcpm <- cpm(x, log=TRUE)
table(rowSums(x$counts==0)==9)
##
## FALSE TRUE
## 22026 5153
keep.exprs <- rowSums(cpm>1)>=3
x <- x[keep.exprs,, keep.lib.sizes=FALSE]
dim(x)
## [1] 14165 9
library(RColorBrewer)
nsamples <- ncol(x)
col <- brewer.pal(nsamples, "Paired")
par(mfrow=c(1,2))
plot(density(lcpm[,1]), col=col[1], lwd=2, ylim=c(0,0.21), las=2,
main="", xlab="")
title(main="A. Raw data", xlab="Log-cpm")
abline(v=0, lty=3)
for (i in 2:nsamples){
den <- density(lcpm[,i])
lines(den$x, den$y, col=col[i], lwd=2)
}
legend("topright", samplenames, text.col=col, bty="n")
lcpm <- cpm(x, log=TRUE)
plot(density(lcpm[,1]), col=col[1], lwd=2, ylim=c(0,0.21), las=2,
main="", xlab="")
title(main="B. Filtered data", xlab="Log-cpm")
abline(v=0, lty=3)
for (i in 2:nsamples){
den <- density(lcpm[,i])
lines(den$x, den$y, col=col[i], lwd=2)
}
legend("topright", samplenames, text.col=col, bty="n")
x <- calcNormFactors(x, method = "TMM")
x$samples$norm.factors
## [1] 0.896 1.035 1.044 1.041 1.032 0.922 0.984 1.083 0.979
x2 <- x
x2$samples$norm.factors <- 1
x2$counts[,1] <- ceiling(x2$counts[,1]*0.05)
x2$counts[,2] <- x2$counts[,2]*5
par(mfrow=c(1,2))
lcpm <- cpm(x2, log=TRUE)
boxplot(lcpm, las=2, col=col, main="")
title(main="A. Example: Unnormalised data",ylab="Log-cpm")
x2 <- calcNormFactors(x2)
x2$samples$norm.factors
## [1] 0.0547 6.1306 1.2293 1.1705 1.2149 1.0562 1.1459 1.2613 1.1170
lcpm <- cpm(x2, log=TRUE)
boxplot(lcpm, las=2, col=col, main="")
title(main="B. Example: Normalised data",ylab="Log-cpm")
lcpm <- cpm(x, log=TRUE)
par(mfrow=c(1,2))
col.group <- group
levels(col.group) <- brewer.pal(nlevels(col.group), "Set1")
col.group <- as.character(col.group)
col.lane <- lane
levels(col.lane) <- brewer.pal(nlevels(col.lane), "Set2")
col.lane <- as.character(col.lane)
plotMDS(lcpm, labels=group, col=col.group)
title(main="A. Sample groups")
plotMDS(lcpm, labels=lane, col=col.lane, dim=c(3,4))
title(main="B. Sequencing lanes")
library(Glimma)
glMDSPlot(lcpm, labels=paste(group, lane, sep="_"),
groups=x$samples[,c(2,5)], launch=FALSE)
Link to Interactive MDS plot
Differential expression analysis
design <- model.matrix(~0+group+lane)
colnames(design) <- gsub("group", "", colnames(design))
design
## Basal LP ML laneL006 laneL008
## 1 0 1 0 0 0
## 2 0 0 1 0 0
## 3 1 0 0 0 0
## 4 1 0 0 1 0
## 5 0 0 1 1 0
## 6 0 1 0 1 0
## 7 1 0 0 1 0
## 8 0 0 1 0 1
## 9 0 1 0 0 1
## attr(,"assign")
## [1] 1 1 1 2 2
## attr(,"contrasts")
## attr(,"contrasts")$group
## [1] "contr.treatment"
##
## attr(,"contrasts")$lane
## [1] "contr.treatment"
contr.matrix <- makeContrasts(
BasalvsLP = Basal-LP,
BasalvsML = Basal - ML,
LPvsML = LP - ML,
levels = colnames(design))
contr.matrix
## Contrasts
## Levels BasalvsLP BasalvsML LPvsML
## Basal 1 1 0
## LP -1 0 1
## ML 0 -1 -1
## laneL006 0 0 0
## laneL008 0 0 0
par(mfrow=c(1,2))
v <- voom(x, design, plot=TRUE)
v
## An object of class "EList"
## $genes
## ENTREZID SYMBOL TXCHROM
## 1 497097 Xkr4 chr1
## 6 27395 Mrpl15 chr1
## 7 18777 Lypla1 chr1
## 9 21399 Tcea1 chr1
## 10 58175 Rgs20 chr1
## 14160 more rows ...
##
## $targets
## files group lib.size norm.factors lane
## 10_6_5_11 GSM1545535_10_6_5_11.txt LP 29409426 0.896 L004
## 9_6_5_11 GSM1545536_9_6_5_11.txt ML 36528591 1.035 L004
## purep53 GSM1545538_purep53.txt Basal 59598629 1.044 L004
## JMS8-2 GSM1545539_JMS8-2.txt Basal 53382070 1.041 L006
## JMS8-3 GSM1545540_JMS8-3.txt ML 78175314 1.032 L006
## JMS8-4 GSM1545541_JMS8-4.txt LP 55762781 0.922 L006
## JMS8-5 GSM1545542_JMS8-5.txt Basal 54115150 0.984 L006
## JMS9-P7c GSM1545544_JMS9-P7c.txt ML 23043111 1.083 L008
## JMS9-P8c GSM1545545_JMS9-P8c.txt LP 19525423 0.979 L008
##
## $E
## Samples
## Tags 10_6_5_11 9_6_5_11 purep53 JMS8-2 JMS8-3 JMS8-4 JMS8-5 JMS9-P7c JMS9-P8c
## 497097 -4.29 -3.87 2.52 3.30 -4.48 -3.99 3.31 -3.20 -5.29
## 27395 3.88 4.40 4.52 4.57 4.32 3.79 3.92 4.35 4.13
## 18777 4.71 5.56 5.40 5.17 5.63 5.08 5.08 5.76 5.15
## 21399 4.78 4.74 5.37 5.13 4.85 4.94 5.16 5.04 4.99
## 58175 3.94 3.29 -1.77 -1.88 2.99 3.36 -2.11 3.14 3.52
## 14160 more rows ...
##
## $weights
## [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9]
## [1,] 1.18 1.18 20.53 20.98 1.77 1.22 21.13 1.18 1.18
## [2,] 20.88 26.56 31.60 29.66 32.56 26.75 29.79 21.90 17.15
## [3,] 28.00 33.70 34.85 34.46 35.15 33.55 34.52 31.44 25.23
## [4,] 27.67 29.60 34.90 34.43 34.84 33.16 34.49 26.14 24.50
## [5,] 19.74 18.66 3.18 2.63 24.19 24.01 2.65 13.15 14.35
## 14160 more rows ...
##
## $design
## Basal LP ML laneL006 laneL008
## 1 0 1 0 0 0
## 2 0 0 1 0 0
## 3 1 0 0 0 0
## 4 1 0 0 1 0
## 5 0 0 1 1 0
## 6 0 1 0 1 0
## 7 1 0 0 1 0
## 8 0 0 1 0 1
## 9 0 1 0 0 1
## attr(,"assign")
## [1] 1 1 1 2 2
## attr(,"contrasts")
## attr(,"contrasts")$group
## [1] "contr.treatment"
##
## attr(,"contrasts")$lane
## [1] "contr.treatment"
vfit <- lmFit(v, design)
vfit <- contrasts.fit(vfit, contrasts=contr.matrix)
efit <- eBayes(vfit)
plotSA(efit)
summary(decideTests(efit))
## BasalvsLP BasalvsML LPvsML
## -1 4127 4338 2895
## 0 5740 5655 8825
## 1 4298 4172 2445
tfit <- treat(vfit, lfc=1)
dt <- decideTests(tfit)
summary(dt)
## BasalvsLP BasalvsML LPvsML
## -1 1417 1512 203
## 0 11030 10895 13780
## 1 1718 1758 182
de.common <- which(dt[,1]!=0 & dt[,2]!=0)
length(de.common)
## [1] 2409
head(tfit$genes$SYMBOL[de.common], n=20)
## [1] "Xkr4" "Rgs20" "Cpa6" "Sulf1" "Eya1"
## [6] "Msc" "Sbspon" "Pi15" "Crispld1" "Kcnq5"
## [11] "Ptpn18" "Arhgef4" "2010300C02Rik" "Aff3" "Npas2"
## [16] "Tbc1d8" "Creg2" "Il1r1" "Il18r1" "Il18rap"
vennDiagram(dt[,1:2], circle.col=c("turquoise", "salmon"))
write.fit(tfit, dt, file="results.txt")
basal.vs.lp <- topTreat(tfit, coef=1, n=Inf)
basal.vs.ml <- topTreat(tfit, coef=2, n=Inf)
head(basal.vs.lp)
## ENTREZID SYMBOL TXCHROM logFC AveExpr t P.Value adj.P.Val
## 12759 12759 Clu chr14 -5.44 8.86 -33.4 3.99e-10 2.7e-06
## 53624 53624 Cldn7 chr11 -5.51 6.30 -32.9 4.50e-10 2.7e-06
## 242505 242505 Rasef chr4 -5.92 5.12 -31.8 6.06e-10 2.7e-06
## 67451 67451 Pkp2 chr16 -5.72 4.42 -30.7 8.01e-10 2.7e-06
## 228543 228543 Rhov chr2 -6.25 5.49 -29.5 1.11e-09 2.7e-06
## 70350 70350 Basp1 chr15 -6.07 5.25 -28.6 1.38e-09 2.7e-06
head(basal.vs.ml)
## ENTREZID SYMBOL TXCHROM logFC AveExpr t P.Value adj.P.Val
## 242505 242505 Rasef chr4 -6.51 5.12 -35.5 2.57e-10 1.92e-06
## 53624 53624 Cldn7 chr11 -5.47 6.30 -32.5 4.98e-10 1.92e-06
## 12521 12521 Cd82 chr2 -4.67 7.07 -31.8 5.80e-10 1.92e-06
## 71740 71740 Pvrl4 chr1 -5.56 5.17 -31.3 6.76e-10 1.92e-06
## 20661 20661 Sort1 chr3 -4.91 6.71 -31.2 6.76e-10 1.92e-06
## 15375 15375 Foxa1 chr12 -5.75 5.63 -28.3 1.49e-09 2.28e-06
plotMD(tfit, column=1, status=dt[,1], main=colnames(tfit)[1],
xlim=c(-8,13))
glMDPlot(tfit, coef=1, status=dt[,1], main=colnames(tfit)[1],
counts=x$counts, samples=colnames(x), anno=x$genes, groups=group,
id.column="ENTREZID", display.columns=c("SYMBOL", "ENTREZID"),
search.by="SYMBOL", launch=FALSE)
Link to Interactive MD plot
library(gplots)
basal.vs.lp.topgenes <- basal.vs.lp$ENTREZID[1:100]
i <- which(v$genes$ENTREZID %in% basal.vs.lp.topgenes)
mycol <- colorpanel(1000,"blue","white","red")
heatmap.2(v$E[i,], scale="row",
labRow=v$genes$SYMBOL[i], labCol=group,
col=mycol, trace="none", density.info="none",
margin=c(8,6), lhei=c(2,10), dendrogram="column")
Gene set testing with camera
load(url("http://bioinf.wehi.edu.au/software/MSigDB/mouse_c2_v5p1.rdata"))
idx <- ids2indices(Mm.c2,id=rownames(v))
cam.BasalvsLP <- camera(v,idx,design,contrast=contr.matrix[,1],
inter.gene.cor=0.01)
head(cam.BasalvsLP,5)
## NGenes Correlation Direction PValue
## LIM_MAMMARY_STEM_CELL_UP 739 0.01 Up 1.13e-18
## LIM_MAMMARY_STEM_CELL_DN 630 0.01 Down 1.57e-15
## ROSTY_CERVICAL_CANCER_PROLIFERATION_CLUSTER 163 0.01 Up 1.44e-13
## SOTIRIOU_BREAST_CANCER_GRADE_1_VS_3_UP 183 0.01 Up 2.18e-13
## LIM_MAMMARY_LUMINAL_PROGENITOR_UP 87 0.01 Down 6.73e-13
## FDR
## LIM_MAMMARY_STEM_CELL_UP 5.36e-15
## LIM_MAMMARY_STEM_CELL_DN 3.71e-12
## ROSTY_CERVICAL_CANCER_PROLIFERATION_CLUSTER 2.26e-10
## SOTIRIOU_BREAST_CANCER_GRADE_1_VS_3_UP 2.58e-10
## LIM_MAMMARY_LUMINAL_PROGENITOR_UP 6.36e-10
cam.BasalvsML <- camera(v,idx,design,contrast=contr.matrix[,2],
inter.gene.cor=0.01)
head(cam.BasalvsML,5)
## NGenes Correlation Direction PValue
## LIM_MAMMARY_STEM_CELL_UP 739 0.01 Up 5.09e-23
## LIM_MAMMARY_STEM_CELL_DN 630 0.01 Down 5.13e-19
## LIM_MAMMARY_LUMINAL_MATURE_DN 166 0.01 Up 8.88e-16
## LIM_MAMMARY_LUMINAL_MATURE_UP 180 0.01 Down 6.29e-13
## ROSTY_CERVICAL_CANCER_PROLIFERATION_CLUSTER 163 0.01 Up 1.68e-12
## FDR
## LIM_MAMMARY_STEM_CELL_UP 2.40e-19
## LIM_MAMMARY_STEM_CELL_DN 1.21e-15
## LIM_MAMMARY_LUMINAL_MATURE_DN 1.40e-12
## LIM_MAMMARY_LUMINAL_MATURE_UP 7.43e-10
## ROSTY_CERVICAL_CANCER_PROLIFERATION_CLUSTER 1.59e-09
cam.LPvsML <- camera(v,idx,design,contrast=contr.matrix[,3],
inter.gene.cor=0.01)
head(cam.LPvsML,5)
## NGenes Correlation Direction PValue FDR
## LIM_MAMMARY_LUMINAL_MATURE_UP 180 0.01 Down 8.50e-14 3.40e-10
## LIM_MAMMARY_LUMINAL_MATURE_DN 166 0.01 Up 1.44e-13 3.40e-10
## LIM_MAMMARY_LUMINAL_PROGENITOR_UP 87 0.01 Up 3.84e-11 6.05e-08
## REACTOME_RESPIRATORY_ELECTRON_TRANSPORT 91 0.01 Down 2.66e-08 3.14e-05
## NABA_CORE_MATRISOME 222 0.01 Up 4.43e-08 4.19e-05
barcodeplot(efit$t[,3], index=idx$LIM_MAMMARY_LUMINAL_MATURE_UP,
index2=idx$LIM_MAMMARY_LUMINAL_MATURE_DN, main="LPvsML")
Software and code used
sessionInfo()
## R version 3.3.0 (2016-05-03)
## Platform: x86_64-apple-darwin13.4.0 (64-bit)
## Running under: OS X 10.9.5 (Mavericks)
##
## locale:
## [1] en_AU.UTF-8/en_AU.UTF-8/en_AU.UTF-8/C/en_AU.UTF-8/en_AU.UTF-8
##
## attached base packages:
## [1] parallel stats4 stats graphics grDevices utils datasets methods
## [9] base
##
## other attached packages:
## [1] gplots_3.0.1 Glimma_1.1.1
## [3] RColorBrewer_1.1-2 Mus.musculus_1.3.1
## [5] TxDb.Mmusculus.UCSC.mm10.knownGene_3.2.2 org.Mm.eg.db_3.3.0
## [7] GO.db_3.3.0 OrganismDbi_1.14.1
## [9] GenomicFeatures_1.24.2 GenomicRanges_1.24.0
## [11] GenomeInfoDb_1.8.1 AnnotationDbi_1.34.3
## [13] IRanges_2.6.0 S4Vectors_0.10.1
## [15] Biobase_2.32.0 BiocGenerics_0.18.0
## [17] edgeR_3.15.0 limma_3.29.5
## [19] knitr_1.13 BiocStyle_2.0.2
## [21] BiocInstaller_1.23.4
##
## loaded via a namespace (and not attached):
## [1] splines_3.3.0 R.utils_2.3.0 gtools_3.5.0
## [4] Formula_1.2-1 latticeExtra_0.6-28 RBGL_1.48.1
## [7] Rsamtools_1.24.0 yaml_2.1.13 RSQLite_1.0.0
## [10] lattice_0.20-33 chron_2.3-47 digest_0.6.9
## [13] XVector_0.12.0 colorspace_1.2-6 htmltools_0.3.5
## [16] Matrix_1.2-6 R.oo_1.20.0 plyr_1.8.3
## [19] DESeq2_1.13.3 XML_3.98-1.4 biomaRt_2.28.0
## [22] genefilter_1.54.2 zlibbioc_1.18.0 xtable_1.8-2
## [25] scales_0.4.0 gdata_2.17.0 BiocParallel_1.6.2
## [28] annotate_1.50.0 ggplot2_2.1.0 SummarizedExperiment_1.2.2
## [31] nnet_7.3-12 survival_2.39-4 magrittr_1.5
## [34] evaluate_0.9 R.methodsS3_1.7.1 foreign_0.8-66
## [37] graph_1.50.0 tools_3.3.0 data.table_1.9.6
## [40] formatR_1.4 stringr_1.0.0 munsell_0.4.3
## [43] locfit_1.5-9.1 cluster_2.0.4 Biostrings_2.40.1
## [46] caTools_1.17.1 grid_3.3.0 RCurl_1.95-4.8
## [49] bitops_1.0-6 rmarkdown_0.9.6 gtable_0.2.0
## [52] DBI_0.4-1 GenomicAlignments_1.8.0 gridExtra_2.2.1
## [55] rtracklayer_1.32.0 Hmisc_3.17-4 KernSmooth_2.23-15
## [58] stringi_1.1.1 Rcpp_0.12.5 geneplotter_1.50.0
## [61] rpart_4.1-10 acepack_1.3-3.3