Versions Compared

Version Old Version 72 New Version Current
Changes made by

Jeffrey E Barrick

Jeffrey E Barrick

Saved on

Key

  • This line was added.
  • This line was removed.
  • Formatting was changed.
Comment: Migrated to Confluence 5.3

...

We've done this several times before, so you should be able to come up with the full command lines if you refer back to the original lesson.

Warning

Be careful we are now mapping single-end reads, so you may have to look at the bowtie help to figure out how to do that!

...

Expand
Please take me through all of the steps...
Please take me through all of the steps...
Code Block
module load bowtie
bowtie-build NC_017544.1.fasta NC_017544.1

Now create a bowtie_commands file that looks like this using nano or another text editor:

Code Block
bowtie -p 3 -S NC_017544.1 SRR034450.fastq -S SRR034450.sam
bowtie -p 3 -S NC_017544.1 SRR034451.fastq -S SRR034451.sam
bowtie -p 3 -S NC_017544.1 SRR034452.fastq -S SRR034452.sam
bowtie -p 3 -S NC_017544.1 SRR034453.fastq -S SRR034453.sam

Remember that there are 12 processors per node on Lonestar, so we choose to use 3 for each of the 4 jobs with the -p 3 option.

Create the launcher script and run it:

Warning

Remember that you cannot qsub from within an idev shell!

Code Block
module load python
launcher_creator.py -n bowtie -q development -cj bowtie_commands -t 0:30:00
qsub bowtie.sge

...

Code Block
titleUsing edgeR
login1$ R
...
> library("edgeR")
> counts = read.delim("gene_counts.tab", header=F, row.names=1)
> colnames(counts) = c("wt1", "mut1", "wt2", "mut2")
> head(counts)
> group <- factor(c("wt","mut","wt","mut"))
> dge = DGEList(counts=counts,group=group)
> dge <- estimateCommonDisp(dge)
> dge <- estimateTagwiseDisp(dge)
> et <- exactTest(dge)
> etp <- topTags(et, n=100000)
> etp$table$logFC = -etp$table$logFC

> pdf("edgeR-MA-plot.pdf")
> plot(
  etp$table$logCPM,
  etp$table$logFC,
  xlim=c(-3, 20), ylim=c(-12, 12), pch=20, cex=.3,
  col = ifelse( etp$table$FDR < .1, "red", "black" ) )
> dev.off()

> write.csv(etp$table, "edgeR-wt-vs-mut.csv")
> q()
Save workspace image? [y/n/c]: n
login1$ head edgeR-wt-vs-mut.csv

Note that the "FC" fold change calculated is initially the reverse of that for the DESeq example for the output here. It is wt relative to mut. To fix this, we put a negative in there for the log fold change.

...