# data used in this workshop can be found at - http://bioinf.wehi.edu.au/RNAseqCaseStudy/
# the data size is 283MB

library(Rsubread)
library(limma)
library(edgeR)
options(digits=2)
targets <- readTargets()
celltype <- factor(targets$CellType)
design <- model.matrix(~celltype)
targets
buildindex(basename="chr1",reference="hg19_chr1.fa")
for(i in 1:4)
  align(index="chr1",readfile1=targets$InputFile[i],output_file=targets$OutputFile[i],unique=TRUE,indels=5)
fcounts <- featureCounts(files=targets$OutputFile,genome="hg")
isexpr <- rowSums(cpm(fcounts$counts) > 10) >= 2
x <- fcounts$counts[isexpr,]
y <- voom(x,design,plot=TRUE)
plotMDS(y,xlim=c(-2.5,2.5))
fit <- eBayes(lmFit(y,design))
topTable(fit,coef=2)