Contents



All big data files (processed read data, assemblies, ...) as well as execution commands and Blast results are also available at the Open Science Framework under accession doi.org/10.17605/OSF.IO/5ZDX4.

S1: Data sets and preprocessing

back to top

Overview of all data sets used for assembly. All reads were quality checked with FASTQC and Prinseq prior assembly. Shown are read statistics before and after trimming. If a data set was sequenced by preserving the strand specificity (ss), this is indicated by the direction of the sequenced read (pairs): F -- forward, R -- reverse.

IDSpeciesDomainEncodingTypeStranded# Reads (raw)Read length (raw)%GC (raw)FastQC (raw)# Reads (trimmed)Read length (trimmed)%GC (trimmed)FastQC (trimmed)
ECOEscherichia coliBacteria Illumina 1.9single-endss F79431309440%SINGLE664370925-9450%SINGLE
CAL3Candida albicansFungi Illumina 1.9paired-endnot115769325138%FORWARD | REVERSE1116852925-5138%FORWARD | REVERSE
ATHArabidopsis thalianaPlant Illumina 1.9single-endnot1691177430-10146%SINGLE1682803625-10146%SINGLE
MMUMus musculusMammal Illumina 1.9paired-endss FR526452387648%FORWARD | REVERSE4332153725-7648%FORWARD | REVERSE
HSAHomo sapiensMammal Illumina 1.9paired-endss FR9754805210153%FORWARD | REVERSE9609311625-10152%FORWARD | REVERSE
EBOV_HSA_3HHomo sapiens + EBOV 3hMammal+Virus Illumina 1.9paired-endnot1720378520-10048%FORWARD | REVERSE1564976325-10047%FORWARD | REVERSE
EBOV_HSA_7HHomo sapiens + EBOV 7hMammal+Virus Illumina 1.9paired-endnot2468852320-10047%FORWARD | REVERSE2226860125-10047%FORWARD | REVERSE
EBOV_HSA_23HHomo sapiens + EBOV 23hMammal+Virus Illumina 1.9paired-endnot2647015520-10046%FORWARD | REVERSE2393025025-10046%FORWARD | REVERSE
HSA_FLUXHomo sapiens simulatedMammal-simulated Illumina 1.9paired-endnot600021783-10050%FORWARD | REVERSE5696218125-10049%FORWARD | REVERSE



S2: Assembly tools

back to top

Overview about the evaluated assembly tools. We obtained the most recent versions in November 2018. We chose Oases and Trans-ABySS as two transcriptome assemblers build on top of the genome assemblers Velvet and ABySS, respectively. To achieve an assembly with Oases or Trans-ABySS one first has to run the underlying genome assembly tool with a range of different k-mers. SOAPdenovo-Trans, build on the principles of SOAPdenovo2, can be used as a stand alone tool for de novo transcriptome assembly based on single k-mer values. Trinity was also especially designed for transcriptome assembly, using only a fixed k-mer value of 25. All these four transcriptome assemblers are designed in the scope of working with RNA-Seq data and are based on de Bruijn graph algorithms. In contrast the Mira assembler works on overlap- graphs and therefore uses no k-mer approach to decompose reads before assembly. It can be run in EST mode to deal with the special requirements of RNA-Seq data. Next to this transcriptome assemblers we also chose one de novo genome assembly tool, SPAdes, originally developed for smaller bacterial-size genomes and single-cell data and also based on de Bruijn graph and multiple k-mer values. MK – Yes, if the tool has an build in multiple kmer approach and automatically merges the output of different kmer runs. aOases was used on top of the de novo genome assembler Velvet (v1.2.10). bSPAdes, originally designed as a de novo genome assembler for single-cell data, was used in RNA-Seq modus (-rna) and single-cell modus (-sc), respectively. cWhen running SPAdes in RNA-Seq modus, two k-mer values are used.

AssemblerVersionMultiple k-mer modePMID
Trinity2.8.4no21572440
Trans-ABySS2.0.1yes20935650
SOAPdenovo-Trans1.03no24532719
Oasesa0.2.08yes22368243
IDBA-Tran1.1.1yes23813001
Bridgerv2014-12-01no25723335
BinPacker1.0no26894997
Shannon0.0.2nobioRxiv
SPAdes-scb3.13.0yes22506599
SPAdes-rnab3.13.0yescbioRxiv



S3: Executed assembly commands

back to top

In the following we provide all executed commands that were used to run the different assemblies. Each assembler listed in Tab. S2 was run on each data set listed in Tab. S1. For details like used k-mers and strand-specificity download the readme files or expand the view for the corresponding data set. The command files still include the execution command of MIRA, however the output was not used in the comparison because MIRA performed worse on all data sets (very short contigs, almost no homology to known transcripts, bad re-mapping rate, ... see manuscript).

Assembly of Escherichia coli (download script)
#######################################
## Eco1 1x 94bp, strand specific F
#######################################

##PARAMETERS SHARED
PROJECT=eco
R1=eco1_formated.fastq
THREADS=48
MEM=400
DIR=~/$PROJECT
FINAL_DIR=~/$PROJECT/final

mkdir -p $DIR
mkdir -p $FINAL_DIR

kmer1=25
kmer2=35
kmer3=45
kmer4=55
kmer5=65


###################################################################################
## TRINITY
mkdir -p $DIR/trinity
Trinity --seqType fq --max_memory $MEM"G" --single $R1 --CPU $THREADS --output $DIR/trinity/ --SS_lib_type F
ln -s $DIR/trinity/Trinity.fasta $FINAL_DIR/trinity.fasta
###################################################################################

###################################################################################
## TRANS-ABYSS
mkdir -p $DIR/transabyss

for kmer in $kmer1 $kmer2 $kmer3 $kmer4 $kmer5; do
name=k${kmer}
assemblydir=$DIR/transabyss/${name}
mkdir -p $assemblydir
finalassembly=${assemblydir}/${name}-final.fa
transabyss -k ${kmer} --se $R1 --SS --outdir ${assemblydir} --name ${name} --threads $THREADS
done

mergedassembly=$DIR/transabyss/merged.fa
finalassembly1=$DIR/transabyss/k"$kmer1"/k"$kmer1"-final.fa
finalassembly2=$DIR/transabyss/k"$kmer2"/k"$kmer2"-final.fa
finalassembly3=$DIR/transabyss/k"$kmer3"/k"$kmer3"-final.fa
finalassembly4=$DIR/transabyss/k"$kmer4"/k"$kmer4"-final.fa
finalassembly5=$DIR/transabyss/k"$kmer5"/k"$kmer5"-final.fa
transabyss-merge --threads $THREADS --mink $kmer1 --maxk $kmer5 --SS --prefixes k${kmer1}. k${kmer2}. k${kmer3}. k${kmer4}. k${kmer5}. --out ${mergedassembly} ${finalassembly1} ${finalassembly2} ${finalassembly3} ${finalassembly4} ${finalassembly5}
ln -s $DIR/transabyss/merged.fa $FINAL_DIR/trans-abyss.fa
###################################################################################

###################################################################################
## SOAP-TRANS
mkdir -p $DIR/soap-trans
CONFIG=$DIR/soap-trans/$PROJECT.config
touch $CONFIG
echo '#maximal read length' > $CONFIG
echo "max_rd_len=94" >> $CONFIG
echo "[LIB]" >> $CONFIG
echo "#maximal read length in this lib" >> $CONFIG
echo "rd_len_cutoff=94" >> $CONFIG
#echo "#average insert size" >> $CONFIG
#echo "avg_ins=200" >> $CONFIG
echo "#if sequence needs to be reversed" >> $CONFIG
echo "reverse_seq=0" >> $CONFIG
echo "#in which part(s) the reads are used" >> $CONFIG
echo "asm_flags=3" >> $CONFIG
#echo "q=#{@reads}" >> $CONFIG #if @type == :SE
echo "#fastq paired end files" >> $CONFIG
echo "q="$R1 >> $CONFIG

kmer_min=$kmer1
kmer_max=$kmer5

# run soap with default settings
mkdir -p $DIR/soap-trans/default
SOAPdenovo-Trans-127mer all -s $DIR/soap-trans/$PROJECT.config -o $DIR/soap-trans/default/k23 -p $THREADS
ln -s $DIR/soap-trans/default/k23.scafSeq $FINAL_DIR/soap-trans-default.fasta
###################################################################################

###################################################################################
## RNA-SPADES
mkdir -p $DIR/rnaspades
rnaspades.py --s1 $R1 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/rnaspades/
ln -s $DIR/rnaspades/transcripts.fasta $FINAL_DIR/rna-spades.fasta
###################################################################################

###################################################################################
## SPADES
mkdir -p $DIR/spades
spades.py --sc --s1 $R1 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/spades/
ln -s $DIR/spades/scaffolds.fasta $FINAL_DIR/spades.fasta
###################################################################################

###################################################################################
## OASES
k_min=$kmer1
k_max=$((kmer5 + 1))
k_step=10
type=-short

mkdir -p $DIR/oases
oases_pipeline.py -m $k_min -M $k_max -s $k_step -o $DIR/oases/ -d " $type -fastq $R1 -strand_specific " -c
ln -s $DIR/oases/Merged/transcripts.fa $FINAL_DIR/oases.fasta
###################################################################################

###################################################################################
## IDBA-TRANS
###this tools assume the paired-end reads are in order (->, <-). If your data is in order (<-, ->), please convert it by yourself.
mkdir $DIR/idba-tran
fq2fa --filter $R1 $DIR/idba-tran/r1.fasta
idba_tran -r $DIR/idba-tran/r1.fasta -o $DIR/idba-tran/ --mink $kmer1 --maxk $kmer5 --step 10 --num_threads $THREADS
ln -s $DIR/idba-tran/contig.fa $FINAL_DIR/idba-tran.fa
###################################################################################

###################################################################################
## BRIDGER
OUT=$DIR/bridger/
Bridger.pl --seqType fq --single $R1 --SS_lib_type F --CPU $THREADS -o $OUT
Bridger.pl --seqType fq --single $R1 --SS_lib_type F --CPU $THREADS -o $OUT # start second time, because first run crashs at some point and then finishs
ln -s $OUT/Bridger.fasta $FINAL_DIR/bridger.fasta
###################################################################################

###################################################################################
## BINPACKER
OUT=$DIR/binpacker/
mkdir -p $OUT
cd $OUT
BinPacker -s fq -p single -u $R1 -m F -o $OUT
cd ..
ln -s $OUT/BinPacker_Out_Dir/BinPacker.fa $FINAL_DIR/binpacker.fasta
###################################################################################

###################################################################################
## SHANNON
mkdir $DIR/shannon
python shannon.py -p $THREADS -o $DIR/shannon --single $R1 --ss
ln -s $DIR/shannon/shannon.fasta $FINAL_DIR/shannon.fasta
###################################################################################

###################################################################################
## MIRA
mkdir $DIR/mira
touch $DIR/mira/manifest.config
echo "project = $PROJECT" > $DIR/mira/manifest.config
echo "job = est,denovo,accurate" >> $DIR/mira/manifest.config
echo "parameters = -DI:trt=/mnt/dessertlocal/projects/transcriptome_assembly/tmp -NW:cnfs=no -NW:cmrnl=no -SK:mmhr=1" >> $DIR/mira/manifest.config
echo "# defining illumina paired end reads" >> $DIR/mira/manifest.config
echo "readgroup = se" >> $DIR/mira/manifest.config
echo "data = $R1" >> $DIR/mira/manifest.config
echo "technology = solexa" >> $DIR/mira/manifest.config
#echo "template_size = 50 1000 autorefine" >> $DIR/mira/manifest.config

mira -c $DIR/mira -t $THREADS $DIR/mira/manifest.config
ln -s $DIR/mira/"$PROJECT"_assembly/"$PROJECT"_d_results/"$PROJECT"_out.unpadded.fasta $FINAL_DIR/mira.fasta
###################################################################################


Assembly of Candida albicans (download script)
#######################################
## CAL, 2x 51bp not strand-specific
#######################################

##PARAMETERS SHARED
PROJECT=cal3
R1=cal3_1_formated.fastq
R2=cal3_2_formated.fastq
THREADS=48
MEM=400
DIR=~/$PROJECT
FINAL_DIR=~/$PROJECT/final

mkdir -p $DIR
mkdir -p $FINAL_DIR

kmer1=21
kmer2=27
kmer3=33
kmer4=39

###################################################################################
## TRINITY
mkdir -p $DIR/trinity
Trinity --seqType fq --max_memory $MEM"G" --left $R1 --right $R2 --CPU $THREADS --output $DIR/trinity/
ln -s $DIR/trinity/Trinity.fasta $FINAL_DIR/trinity.fasta
###################################################################################

###################################################################################
## TRANS-ABYSS
mkdir -p $DIR/transabyss

for kmer in $kmer1 $kmer2 $kmer3 $kmer4; do
name=k${kmer}
assemblydir=$DIR/transabyss/${name}
mkdir -p $assemblydir
finalassembly=${assemblydir}/${name}-final.fa
transabyss -k ${kmer} --pe $R1 $R2 --outdir ${assemblydir} --name ${name} --threads $THREADS
done

mergedassembly=$DIR/transabyss/merged.fa
finalassembly1=$DIR/transabyss/k"$kmer1"/k"$kmer1"-final.fa
finalassembly2=$DIR/transabyss/k"$kmer2"/k"$kmer2"-final.fa
finalassembly3=$DIR/transabyss/k"$kmer3"/k"$kmer3"-final.fa
finalassembly4=$DIR/transabyss/k"$kmer4"/k"$kmer4"-final.fa
transabyss-merge --threads $THREADS --mink $kmer1 --maxk $kmer4 --SS --prefixes k${kmer1}. k${kmer2}. k${kmer3}. k${kmer4}. --out ${mergedassembly} ${finalassembly1} ${finalassembly2} ${finalassembly3} ${finalassembly4}
ln -s $DIR/transabyss/merged.fa $FINAL_DIR/trans-abyss.fasta
###################################################################################

###################################################################################
## SOAP-TRANS
mkdir -p $DIR/soap-trans
CONFIG=$DIR/soap-trans/test.config
touch $CONFIG
echo '#maximal read length' > $CONFIG
echo "max_rd_len=51" >> $CONFIG
echo "[LIB]" >> $CONFIG
echo "#maximal read length in this lib" >> $CONFIG
echo "rd_len_cutoff=51" >> $CONFIG
#echo "#average insert size" >> $CONFIG
#echo "avg_ins=200" >> $CONFIG
echo "#if sequence needs to be reversed" >> $CONFIG
echo "reverse_seq=0" >> $CONFIG
echo "#in which part(s) the reads are used" >> $CONFIG
echo "asm_flags=3" >> $CONFIG
echo "#fastq paired end files" >> $CONFIG
echo "q1="$R1 >> $CONFIG
echo "q2="$R2 >> $CONFIG

# run soap with defualt settings
mkdir -p $DIR/soap-trans/default
SOAPdenovo-Trans-127mer all -s $DIR/soap-trans/test.config -o $DIR/soap-trans/default/k23 -p $THREADS
ln -s $DIR/soap-trans/default/k23.scafSeq $FINAL_DIR/soap-trans-default.fasta
###################################################################################

###################################################################################
## RNA-SPADES

##ATTENTION CHANGED KMER BECAUS DEFAULT 55 IS TO BIG
mkdir -p $DIR/rnaspades
rnaspades.py -k 25 --pe1-1 $R1 --pe1-2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/rnaspades/
ln -s $DIR/rnaspades/transcripts.fasta $FINAL_DIR/rna-spades.fasta
###################################################################################

###################################################################################
## SPADES
mkdir -p $DIR/spades
spades.py --sc --pe1-1 $R1 --pe1-2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/spades/
ln -s $DIR/spades/scaffolds.fasta $FINAL_DIR/spades.fasta
##################################################################################

###################################################################################
## OASES
k_min=$kmer1
k_max=$((kmer4 + 1))
k_step=6
type=-shortPaired

mkdir -p $DIR/oases
oases_pipeline.py -m $k_min -M $k_max -s $k_step -o $DIR/oases/ -d " $type -fastq $R1 $R2 " -c
ln -s $DIR/oases/Merged/transcripts.fa $FINAL_DIR/oases.fasta
###################################################################################

###################################################################################
## IDBA-TRANS
###this tools assume the paired-end reads are in order (->, <-). If your data is in order (<-, ->), please convert it by yourself.
mkdir $DIR/idba-tran
fq2fa --merge --filter $R1 $R2 $DIR/idba-tran/r12.fasta
idba_tran -r $DIR/idba-tran/r12.fasta -o $DIR/idba-tran/ --mink $kmer1 --maxk $kmer4 --step 6 --num_threads $THREADS
ln -s $DIR/idba-tran/contig.fa $FINAL_DIR/idba-tran.fasta
###################################################################################

###################################################################################
## BRIDGER
OUT=$DIR/bridger/
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT # start second time, because first run crashs at some point and then finishs
ln -s $OUT/Bridger.fasta $FINAL_DIR/bridger.fasta
###################################################################################

###################################################################################
## BINPACKER
OUT=$DIR/binpacker/
mkdir -p $OUT
BinPacker -s fq -p pair -l $R1 -r $R2 -o $OUT
ln -s $OUT/BinPacker_Out_Dir/BinPacker.fa $FINAL_DIR/binpacker.fasta
###################################################################################

###################################################################################
## SHANNON
mkdir $DIR/shannon
python shannon.py -p $THREADS -o $DIR/shannon --left $R1 --right $R2
ln -s $DIR/shannon/shannon.fasta $FINAL_DIR/shannon.fasta
###################################################################################

###################################################################################
## MIRA
mkdir $DIR/mira
touch $DIR/mira/manifest.config
echo "project = $PROJECT" > $DIR/mira/manifest.config
echo "job = est,denovo,accurate" >> $DIR/mira/manifest.config
echo "parameters = -DI:trt=/mnt/dessertlocal/projects/transcriptome_assembly/tmp -NW:cnfs=no -NW:cmrnl=no -SK:mmhr=1" >> $DIR/mira/manifest.config
echo "# defining illumina paired end reads" >> $DIR/mira/manifest.config
echo "readgroup = pe" >> $DIR/mira/manifest.config
echo "data = $R1 $R2" >> $DIR/mira/manifest.config
echo "technology = solexa" >> $DIR/mira/manifest.config
echo "template_size = 50 1000 autorefine" >> $DIR/mira/manifest.config
echo "segment_placement = ---> <---" >> $DIR/mira/manifest.config

mira -c $DIR/mira -t $THREADS $DIR/mira/manifest.config
ln -s $DIR/mira/"$PROJECT"_assembly/"$PROJECT"_d_results/"$PROJECT"_out.unpadded.fasta $FINAL_DIR/mira.fasta
###################################################################################


Assembly of Arabidopsis thaliana (download script)
#######################################
## Arabidopsis 1x 100bp, not strand specific
#######################################

##PARAMETERS SHARED
PROJECT=ath
R1=ath_formated.fastq
THREADS=48
MEM=400
DIR=~/$PROJECT
FINAL_DIR=~/$PROJECT/final

mkdir -p $DIR
mkdir -p $FINAL_DIR

kmer1=25
kmer2=35
kmer3=45
kmer4=55
kmer5=65

###################################################################################
## TRINITY
mkdir -p $DIR/trinity
Trinity --seqType fq --max_memory $MEM"G" --single $R1 --CPU $THREADS --output $DIR/trinity/
ln -s $DIR/trinity/Trinity.fasta $FINAL_DIR/trinity.fasta
###################################################################################

###################################################################################
## TRANS-ABYSS
mkdir -p $DIR/transabyss

for kmer in $kmer1 $kmer2 $kmer3 $kmer4 $kmer5; do
name=k${kmer}
assemblydir=$DIR/transabyss/${name}
mkdir -p $assemblydir
finalassembly=${assemblydir}/${name}-final.fa
transabyss -k ${kmer} --se $R1 --outdir ${assemblydir} --name ${name} --threads $THREADS
done

mergedassembly=$DIR/transabyss/merged.fa
finalassembly1=$DIR/transabyss/k"$kmer1"/k"$kmer1"-final.fa
finalassembly2=$DIR/transabyss/k"$kmer2"/k"$kmer2"-final.fa
finalassembly3=$DIR/transabyss/k"$kmer3"/k"$kmer3"-final.fa
finalassembly4=$DIR/transabyss/k"$kmer4"/k"$kmer4"-final.fa
finalassembly5=$DIR/transabyss/k"$kmer5"/k"$kmer5"-final.fa
transabyss-merge --threads $THREADS --mink $kmer1 --maxk $kmer5 --prefixes k${kmer1}. k${kmer2}. k${kmer3}. k${kmer4}. k${kmer5}. --out ${mergedassembly} ${finalassembly1} ${finalassembly2} ${finalassembly3} ${finalassembly4} ${finalassembly5}
ln -s $DIR/transabyss/merged.fa $FINAL_DIR/trans-abyss.fasta
###################################################################################

###################################################################################
## SOAP-TRANS
mkdir -p $DIR/soap-trans
CONFIG=$DIR/soap-trans/$PROJECT.config
touch $CONFIG
echo '#maximal read length' > $CONFIG
echo "max_rd_len=100" >> $CONFIG
echo "[LIB]" >> $CONFIG
echo "#maximal read length in this lib" >> $CONFIG
echo "rd_len_cutoff=100" >> $CONFIG
#echo "#average insert size" >> $CONFIG
#echo "avg_ins=200" >> $CONFIG
echo "#if sequence needs to be reversed" >> $CONFIG
echo "reverse_seq=0" >> $CONFIG
echo "#in which part(s) the reads are used" >> $CONFIG
echo "asm_flags=3" >> $CONFIG
#echo "q=#{@reads}" >> $CONFIG #if @type == :SE
echo "#fastq single end files" >> $CONFIG
echo "q="$R1 >> $CONFIG

kmer_min=$kmer1
kmer_max=$kmer5

# run soap with default settings
mkdir -p $DIR/soap-trans/default
SOAPdenovo-Trans-127mer all -s $DIR/soap-trans/$PROJECT.config -o $DIR/soap-trans/default/k23 -p $THREADS
ln -s $DIR/soap-trans/default/k23.scafSeq $FINAL_DIR/soap-trans-default.fasta
###################################################################################

###################################################################################
## RNA-SPADES
mkdir -p $DIR/rnaspades
rnaspades.py --s1 $R1 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/rnaspades/
ln -s $DIR/rnaspades/transcripts.fasta $FINAL_DIR/rna-spades.fasta
###################################################################################

###################################################################################
## SPADES
mkdir -p $DIR/spades
spades.py --sc --s1 $R1 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/spades/
ln -s $DIR/spades/scaffolds.fasta $FINAL_DIR/spades.fasta
###################################################################################

###################################################################################
## OASES
k_min=$kmer1
k_max=$((kmer5 + 1))
k_step=10
type=-short

mkdir -p $DIR/oases
oases_pipeline.py -m $k_min -M $k_max -s $k_step -o $DIR/oases/ -d " $type -fastq $R1 " -c
ln -s $DIR/oases/Merged/transcripts.fa $FINAL_DIR/oases.fasta
###################################################################################

###################################################################################
## IDBA-TRANS
###this tools assume the paired-end reads are in order (->, <-). If your data is in order (<-, ->), please convert it by yourself.
mkdir $DIR/idba-tran
fq2fa --filter $R1 $DIR/idba-tran/r1.fasta
idba_tran -r $DIR/idba-tran/r1.fasta -o $DIR/idba-tran/ --mink $kmer1 --maxk $kmer5 --step 10 --num_threads $THREADS
ln -s $DIR/idba-tran/contig.fa $FINAL_DIR/idba-tran.fasta
###################################################################################

###################################################################################
## BRIDGER
OUT=$DIR/bridger/
Bridger.pl --seqType fq --single $R1 --CPU $THREADS -o $OUT
Bridger.pl --seqType fq --single $R1 --CPU $THREADS -o $OUT # start second time, because first run crashs at some point and then finishs
ln -s $OUT/Bridger.fasta $FINAL_DIR/bridger.fasta
###################################################################################

###################################################################################
## BINPACKER
OUT=$DIR/binpacker/
mkdir -p $OUT
cd $OUT
BinPacker -s fq -p single -u $R1 -o $OUT
cd ..
ln -s $OUT/BinPacker_Out_Dir/BinPacker.fa $FINAL_DIR/binpacker.fasta
###################################################################################

###################################################################################
## SHANNON
mkdir $DIR/shannon
python shannon.py -p $THREADS -o $DIR/shannon --single $R1
ln -s $DIR/shannon/shannon.fasta $FINAL_DIR/shannon.fasta
###################################################################################

###################################################################################
## MIRA
mkdir $DIR/mira
touch $DIR/mira/manifest.config
echo "project = $PROJECT" > $DIR/mira/manifest.config
echo "job = est,denovo,accurate" >> $DIR/mira/manifest.config
echo "parameters = -DI:trt=/mnt/dessertlocal/projects/transcriptome_assembly/tmp -NW:cnfs=no -NW:cmrnl=no -SK:mmhr=1" >> $DIR/mira/manifest.config
echo "# defining illumina paired end reads" >> $DIR/mira/manifest.config
echo "readgroup = se" >> $DIR/mira/manifest.config
echo "data = $R1" >> $DIR/mira/manifest.config
echo "technology = solexa" >> $DIR/mira/manifest.config

mira -c $DIR/mira -t $THREADS $DIR/mira/manifest.config
ln -s $DIR/mira/"$PROJECT"_assembly/"$PROJECT"_d_results/"$PROJECT"_out.unpadded.fasta $FINAL_DIR/mira.fasta
###################################################################################


Assembly of Mus musculus (download script)
#######################################
## MMU TRINITY TEST DATA, 2x 76bp, strand specific RF
#######################################

##PARAMETERS SHARED
PROJECT=mmu
R1=mmu_1_formated.fastq
R2=mmu_2_formated.fastq
THREADS=48
MEM=400
DIR=~/$PROJECT
FINAL_DIR=~/$PROJECT/final

mkdir -p $DIR
mkdir -p $FINAL_DIR

kmer1=25
kmer2=35
kmer3=45
kmer4=55

###################################################################################
## TRINITY
mkdir -p $DIR/trinity
Trinity --seqType fq --max_memory $MEM"G" --left $R1 --right $R2 --CPU $THREADS --output $DIR/trinity/ --SS_lib_type FR
ln -s $DIR/trinity/Trinity.fasta $FINAL_DIR/trinity.fasta
###################################################################################

###################################################################################
## TRANS-ABYSS
mkdir -p $DIR/transabyss

for kmer in $kmer1 $kmer2 $kmer3 $kmer4; do
name=k${kmer}
assemblydir=$DIR/transabyss/${name}
mkdir -p $assemblydir
finalassembly=${assemblydir}/${name}-final.fa
transabyss -k ${kmer} --pe $R1 $R2 --SS --outdir ${assemblydir} --name ${name} --threads $THREADS
done

mergedassembly=$DIR/transabyss/merged.fa
finalassembly1=$DIR/transabyss/k"$kmer1"/k"$kmer1"-final.fa
finalassembly2=$DIR/transabyss/k"$kmer2"/k"$kmer2"-final.fa
finalassembly3=$DIR/transabyss/k"$kmer3"/k"$kmer3"-final.fa
finalassembly4=$DIR/transabyss/k"$kmer4"/k"$kmer4"-final.fa
transabyss-merge --threads $THREADS --mink $kmer1 --maxk $kmer4 --SS --prefixes k${kmer1}. k${kmer2}. k${kmer3}. k${kmer4}. --out ${mergedassembly} ${finalassembly1} ${finalassembly2} ${finalassembly3} ${finalassembly4}
ln -s $DIR/transabyss/merged.fa $FINAL_DIR/trans-abyss.fasta
###################################################################################

###################################################################################
## SOAP-TRANS
mkdir -p $DIR/soap-trans
CONFIG=$DIR/soap-trans/test.config
touch $CONFIG
echo '#maximal read length' > $CONFIG
echo "max_rd_len=76" >> $CONFIG
echo "[LIB]" >> $CONFIG
echo "#maximal read length in this lib" >> $CONFIG
echo "rd_len_cutoff=76" >> $CONFIG
#echo "#average insert size" >> $CONFIG
#echo "avg_ins=200" >> $CONFIG
echo "#if sequence needs to be reversed" >> $CONFIG
echo "reverse_seq=1" >> $CONFIG
echo "#in which part(s) the reads are used" >> $CONFIG
echo "asm_flags=3" >> $CONFIG
#echo "q=#{@reads}" >> $CONFIG #if @type == :SE
echo "#fastq paired end files" >> $CONFIG
echo "q1="$R1 >> $CONFIG
echo "q2="$R2 >> $CONFIG

kmer_min=$kmer1
kmer_max=$kmer4

# run soap with defualt settings
mkdir -p $DIR/soap-trans/default
SOAPdenovo-Trans-127mer all -s $DIR/soap-trans/test.config -o $DIR/soap-trans/default/k23 -p $THREADS
ln -s $DIR/soap-trans/default/k23.scafSeq $FINAL_DIR/soap-trans-default.fasta
###################################################################################

###################################################################################
## RNA-SPADES
mkdir -p $DIR/rnaspades
rnaspades.py --pe1-1 $R1 --pe1-2 $R2 --pe1-fr -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/rnaspades/
ln -s $DIR/rnaspades/transcripts.fasta $FINAL_DIR/rna-spades.fasta
###################################################################################

###################################################################################
## SPADES
mkdir -p $DIR/spades
spades.py --sc --pe1-1 $R1 --pe1-2 $R2 --pe1-fr -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/spades/
ln -s $DIR/spades/scaffolds.fasta $FINAL_DIR/spades.fasta
##################################################################################

###################################################################################
## OASES
k_min=$kmer1
k_max=$((kmer4 + 1))
k_step=10
type=-shortPaired

mkdir -p $DIR/oases
oases_pipeline.py -m $k_min -M $k_max -s $k_step -o $DIR/oases/ -d " $type -fastq $R1 $R2 -strand_specific -separate " -c
ln -s $DIR/oases/Merged/transcripts.fa $FINAL_DIR/oases.fasta
###################################################################################

###################################################################################
## IDBA-TRANS
###this tools assume the paired-end reads are in order (->, <-). If your data is in order (<-, ->), please convert it by yourself.
mkdir $DIR/idba-tran
R1=/mnt/dessertlocal/projects/transcriptome_assembly/data/prinseq/mmu_1_formated_reversed.fastq
R2=/mnt/dessertlocal/projects/transcriptome_assembly/data/prinseq/mmu_2_formated_reversed.fastq
fq2fa --merge --filter $R1 $R2 $DIR/idba-tran/r12.fasta
idba_tran -r $DIR/idba-tran/r12.fasta -o $DIR/idba-tran/ --mink $kmer1 --maxk $kmer4 --step 10 --num_threads $THREADS
ln -s $DIR/idba-tran/contig.fa $FINAL_DIR/idba-tran.fasta
R1=/mnt/dessertlocal/projects/transcriptome_assembly/data/prinseq/mmu_1_formated.fastq
R2=/mnt/dessertlocal/projects/transcriptome_assembly/data/prinseq/mmu_2_formated.fastq
###################################################################################

###################################################################################
## BRIDGER
#Error! try unstranded mode!
OUT=$DIR/bridger/
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT # start second time, because first run crashs at some point and then finishs
ln -s $OUT/Bridger.fasta $FINAL_DIR/bridger.fasta
###################################################################################

###################################################################################
## BINPACKER
OUT=$DIR/binpacker/
mkdir -p $OUT
BinPacker -s fq -p pair -l $R1 -r $R2 -m RF -o $OUT
ln -s $OUT/BinPacker_Out_Dir/BinPacker.fa $FINAL_DIR/binpacker.fasta
###################################################################################

###################################################################################
## SHANNON
mkdir $DIR/shannon
python shannon.py -p $THREADS -o $DIR/shannon --left $R1 --right $R2 --ss
ln -s $DIR/shannon/shannon.fasta $FINAL_DIR/shannon.fasta
###################################################################################

###################################################################################
## MIRA
mkdir $DIR/mira
touch $DIR/mira/manifest.config
echo "project = $PROJECT" > $DIR/mira/manifest.config
echo "job = est,denovo,accurate" >> $DIR/mira/manifest.config
echo "parameters = -DI:trt=/mnt/dessertlocal/projects/transcriptome_assembly/tmp -NW:cnfs=no -NW:cmrnl=no -SK:mmhr=1" >> $DIR/mira/manifest.config
# if @type == :SE
# config << "# defining illumina single end reads"
# config << "readgroup = se"
# config << "data = #{@reads}"
# config << "technology = #{@technology}"
# end
echo "# defining illumina paired end reads" >> $DIR/mira/manifest.config
echo "readgroup = pe" >> $DIR/mira/manifest.config
echo "data = $R1 $R2" >> $DIR/mira/manifest.config
echo "technology = solexa" >> $DIR/mira/manifest.config
echo "template_size = 50 1000 autorefine" >> $DIR/mira/manifest.config
#echo "segment_placement = ---> <---" >> $DIR/mira/manifest.config
echo "segment_placement = <--- --->" >> $DIR/mira/manifest.config

mira -c $DIR/mira -t $THREADS $DIR/mira/manifest.config
ln -s $DIR/mira/"$PROJECT"_assembly/"$PROJECT"_d_results/"$PROJECT"_out.unpadded.fasta $FINAL_DIR/mira.fasta
###################################################################################


Assembly of Homo sapiens (download script)
#######################################
## HSA SRA, 2x 100bp, strand specific RF
#######################################

##PARAMETERS SHARED
PROJECT=hsa
R1=hsa_1.fastq
R2=hsa_2.fastq
THREADS=48
MEM=400
DIR=~/$PROJECT
FINAL_DIR=~/$PROJECT/final

mkdir -p $DIR
mkdir -p $FINAL_DIR

kmer1=25
kmer2=35
kmer3=45
kmer4=55
kmer5=65

###################################################################################
## TRINITY
mkdir -p $DIR/trinity
Trinity --seqType fq --max_memory $MEM"G" --left $R1 --right $R2 --CPU $THREADS --output $DIR/trinity/ --SS_lib_type FR
ln -s $DIR/trinity/Trinity.fasta $FINAL_DIR/trinity.fasta
###################################################################################

###################################################################################
## TRANS-ABYSS
mkdir -p $DIR/transabyss

for kmer in $kmer1 $kmer2 $kmer3 $kmer4 $kmer5; do
name=k${kmer}
assemblydir=$DIR/transabyss/${name}
mkdir -p $assemblydir
finalassembly=${assemblydir}/${name}-final.fa
transabyss -k ${kmer} --pe $R1 $R2 --SS --outdir ${assemblydir} --name ${name} --threads $THREADS
done

mergedassembly=$DIR/transabyss/merged.fa
finalassembly1=$DIR/transabyss/k"$kmer1"/k"$kmer1"-final.fa
finalassembly2=$DIR/transabyss/k"$kmer2"/k"$kmer2"-final.fa
finalassembly3=$DIR/transabyss/k"$kmer3"/k"$kmer3"-final.fa
finalassembly4=$DIR/transabyss/k"$kmer4"/k"$kmer4"-final.fa
finalassembly5=$DIR/transabyss/k"$kmer5"/k"$kmer5"-final.fa
transabyss-merge --threads $THREADS --mink $kmer1 --maxk $kmer5 --SS --prefixes k${kmer1}. k${kmer2}. k${kmer3}. k${kmer4}. k${kmer5}. --out ${mergedassembly} ${finalassembly1} ${finalassembly2} ${finalassembly3} ${finalassembly4} ${finalassembly5}
ln -s $DIR/transabyss/merged.fa $FINAL_DIR/trans-abyss.fasta
###################################################################################

###################################################################################
## SOAP-TRANS
mkdir -p $DIR/soap-trans
CONFIG=$DIR/soap-trans/test.config
touch $CONFIG
echo '#maximal read length' > $CONFIG
echo "max_rd_len=100" >> $CONFIG
echo "[LIB]" >> $CONFIG
echo "#maximal read length in this lib" >> $CONFIG
echo "rd_len_cutoff=100" >> $CONFIG
#echo "#average insert size" >> $CONFIG
#echo "avg_ins=200" >> $CONFIG
echo "#if sequence needs to be reversed" >> $CONFIG
echo "reverse_seq=1" >> $CONFIG
echo "#in which part(s) the reads are used" >> $CONFIG
echo "asm_flags=3" >> $CONFIG
#echo "q=#{@reads}" >> $CONFIG #if @type == :SE
echo "#fastq paired end files" >> $CONFIG
echo "q1="$R1 >> $CONFIG
echo "q2="$R2 >> $CONFIG

kmer_min=$kmer1
kmer_max=$kmer5

# run soap with defualt settings
mkdir -p $DIR/soap-trans/default
SOAPdenovo-Trans-127mer all -s $DIR/soap-trans/test.config -o $DIR/soap-trans/default/k23 -p $THREADS
ln -s $DIR/soap-trans/default/k23.scafSeq $FINAL_DIR/soap-trans-default.fasta
###################################################################################

###################################################################################
## RNA-SPADES
mkdir -p $DIR/rnaspades
rnaspades.py --pe1-1 $R1 --pe1-2 $R2 --pe1-fr -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/rnaspades/
ln -s $DIR/rnaspades/transcripts.fasta $FINAL_DIR/rna-spades.fasta
###################################################################################

###################################################################################
## SPADES
mkdir -p $DIR/spades
spades.py --sc --pe1-1 $R1 --pe1-2 $R2 --pe1-fr -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/spades/
ln -s $DIR/spades/scaffolds.fasta $FINAL_DIR/spades.fasta
###################################################################################

###################################################################################
## OASES
k_min=$kmer1
k_max=$((kmer5 + 1))
k_step=10
type=-shortPaired

mkdir -p $DIR/oases
oases_pipeline.py -m $k_min -M $k_max -s $k_step -o $DIR/oases/ -d " $type -fastq $R1 $R2 -strand_specific -separate " -c
ln -s $DIR/oases/Merged/transcripts.fa $FINAL_DIR/oases.fasta
###################################################################################

###################################################################################
## IDBA-TRANS
###this tools assume the paired-end reads are in order (->, <-). If your data is in order (<-, ->), please convert it by yourself.
mkdir $DIR/idba-tran
R1=/mnt/dessertlocal/projects/transcriptome_assembly/data/prinseq/hsa_1_reversed.fastq
R2=/mnt/dessertlocal/projects/transcriptome_assembly/data/prinseq/hsa_2_reversed.fastq
fq2fa --merge --filter $R1 $R2 $DIR/idba-tran/r12.fasta
idba_tran -r $DIR/idba-tran/r12.fasta -o $DIR/idba-tran/ --mink $kmer1 --maxk $kmer5 --step 10 --num_threads $THREADS
ln -s $DIR/idba-tran/contig.fa $FINAL_DIR/idba-tran.fasta
R1=/mnt/dessertlocal/projects/transcriptome_assembly/data/prinseq/hsa_1.fastq
R2=/mnt/dessertlocal/projects/transcriptome_assembly/data/prinseq/hsa_2.fastq
###################################################################################

###################################################################################
## BRIDGER
# Error, try unstranded
OUT=$DIR/bridger/
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT # start second time, because first run crashs at some point and then finishs
ln -s $OUT/Bridger.fasta $FINAL_DIR/bridger.fasta
###################################################################################

###################################################################################
## BINPACKER
OUT=$DIR/binpacker/
mkdir -p $OUT
BinPacker -s fq -p pair -l $R1 -r $R2 -m RF -o $OUT
ln -s $OUT/BinPacker_Out_Dir/BinPacker.fa $FINAL_DIR/binpacker.fasta
###################################################################################

###################################################################################
## SHANNON
mkdir $DIR/shannon
python /mnt/prostlocal/programs/shannon/0.0.2/shannon.py -p $THREADS -o $DIR/shannon --left $R1 --right $R2 --ss
ln -s $DIR/shannon/shannon.fasta $FINAL_DIR/shannon.fasta
###################################################################################

###################################################################################
## MIRA
mkdir $DIR/mira
touch $DIR/mira/manifest.config
echo "project = $PROJECT" > $DIR/mira/manifest.config
echo "job = est,denovo,accurate" >> $DIR/mira/manifest.config
echo "parameters = -DI:trt=/mnt/dessertlocal/projects/transcriptome_assembly/tmp -NW:cnfs=no -NW:cmrnl=no -SK:mmhr=1" >> $DIR/mira/manifest.config
# if @type == :SE
# config << "# defining illumina single end reads"
# config << "readgroup = se"
# config << "data = #{@reads}"
# config << "technology = #{@technology}"
# end
echo "# defining illumina paired end reads" >> $DIR/mira/manifest.config
echo "readgroup = pe" >> $DIR/mira/manifest.config
echo "data = $R1 $R2" >> $DIR/mira/manifest.config
echo "technology = solexa" >> $DIR/mira/manifest.config
echo "template_size = 50 1000 autorefine" >> $DIR/mira/manifest.config
#echo "segment_placement = ---> <---" >> $DIR/mira/manifest.config
echo "segment_placement = <--- --->" >> $DIR/mira/manifest.config

mira -c $DIR/mira -t $THREADS $DIR/mira/manifest.config
ln -s $DIR/mira/"$PROJECT"_assembly/"$PROJECT"_d_results/"$PROJECT"_out.unpadded.fasta $FINAL_DIR/mira.fasta
###################################################################################


Assembly of Homo sapiens + EBOV 3h (download script)
#######################################
## EBOV HSA 3h
## paired-end, not strand specific, 100bp
#######################################

##PARAMETERS SHARED
PROJECT=ebov_hsa_3h
R1=ebov_hsa_3h_1.fastq
R2=ebov_hsa_3h_2.fastq
THREADS=48
MEM=400
DIR=~/$PROJECT
FINAL_DIR=~/$PROJECT/final

mkdir -p $DIR
mkdir -p $FINAL_DIR

kmer1=25
kmer2=29
kmer3=33
kmer4=37
kmer5=41


###################################################################################
## TRINITY
mkdir -p $DIR/trinity
Trinity --seqType fq --max_memory $MEM"G" --left $R1 --right $R2 --CPU $THREADS --output $DIR/trinity/
ln -s $DIR/trinity/Trinity.fasta $FINAL_DIR/trinity.fasta
###################################################################################

###################################################################################
## TRANS-ABYSS
mkdir -p $DIR/transabyss

for kmer in $kmer1 $kmer2 $kmer3 $kmer4 $kmer5; do
name=k${kmer}
assemblydir=$DIR/transabyss/${name}
mkdir -p $assemblydir
finalassembly=${assemblydir}/${name}-final.fa
transabyss -k ${kmer} --pe $R1 $R2 --outdir ${assemblydir} --name ${name} --threads $THREADS
done

mergedassembly=$DIR/transabyss/merged.fa
finalassembly1=$DIR/transabyss/k"$kmer1"/k"$kmer1"-final.fa
finalassembly2=$DIR/transabyss/k"$kmer2"/k"$kmer2"-final.fa
finalassembly3=$DIR/transabyss/k"$kmer3"/k"$kmer3"-final.fa
finalassembly4=$DIR/transabyss/k"$kmer4"/k"$kmer4"-final.fa
finalassembly5=$DIR/transabyss/k"$kmer5"/k"$kmer5"-final.fa
transabyss-merge --threads $THREADS --mink 25 --maxk 41 --prefixes k${kmer1}. k${kmer2}. k${kmer3}. k${kmer4}. k${kmer5}. --out ${mergedassembly} ${finalassembly1} ${finalassembly2} ${finalassembly3} ${finalassembly4} ${finalassembly5}
ln -s $DIR/transabyss/merged.fa $FINAL_DIR/trans-abyss.fasta
###################################################################################

###################################################################################
## SOAP-TRANS
mkdir -p $DIR/soap-trans
CONFIG=$DIR/soap-trans/test.config
touch $CONFIG
echo '#maximal read length' > $CONFIG
echo "max_rd_len=100" >> $CONFIG
echo "[LIB]" >> $CONFIG
echo "#maximal read length in this lib" >> $CONFIG
echo "rd_len_cutoff=100" >> $CONFIG
#echo "#average insert size" >> $CONFIG
#echo "avg_ins=200" >> $CONFIG
echo "#if sequence needs to be reversed" >> $CONFIG
echo "reverse_seq=0" >> $CONFIG
echo "#in which part(s) the reads are used" >> $CONFIG
echo "asm_flags=3" >> $CONFIG
#echo "q=#{@reads}" >> $CONFIG #if @type == :SE
echo "#fastq paired end files" >> $CONFIG
echo "q1="$R1 >> $CONFIG
echo "q2="$R2 >> $CONFIG

kmer_min=$kmer1
kmer_max=$kmer5

# run soap with defualt settings
mkdir -p $DIR/soap-trans/default
SOAPdenovo-Trans-127mer all -s $DIR/soap-trans/test.config -o $DIR/soap-trans/default/k23 -p $THREADS
ln -s $DIR/soap-trans/default/k23.scafSeq $FINAL_DIR/soap-trans-default.fasta
###################################################################################

###################################################################################
## RNA-SPADES
mkdir -p $DIR/rnaspades/
rnaspades.py -1 $R1 -2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/rnaspades/
ln -s $DIR/rnaspades/transcripts.fasta $FINAL_DIR/rna-spades.fasta
###################################################################################

###################################################################################
## SPADES
mkdir -p $DIR/spades
spades.py --sc -1 $R1 -2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/spades/
ln -s $DIR/spades/scaffolds.fasta $FINAL_DIR/spades.fasta
###################################################################################

###################################################################################
## OASES
k_min=$kmer1
k_max=$((kmer5 + 1))
k_step=4
type=-shortPaired

mkdir -p $DIR/oases
oases_pipeline.py -m $k_min -M $k_max -s $k_step -o $DIR/oases/ -d " $type -fastq $R1 $R2 -separate " -c
ln -s $DIR/oases/Merged/transcripts.fa $FINAL_DIR/oases.fasta
###################################################################################

###################################################################################
## IDBA-TRANS
###this tools assume the paired-end reads are in order (->, <-). If your data is in order (<-, ->), please convert it by yourself.
mkdir $DIR/idba-tran
fq2fa --merge --filter $R1 $R2 $DIR/idba-tran/r12.fasta
idba_tran -r $DIR/idba-tran/r12.fasta -o $DIR/idba-tran/ --mink 25 --maxk 41 --step 4
ln -s $DIR/idba-tran/contig.fa $FINAL_DIR/idba-tran.fasta
###################################################################################

###################################################################################
## BINPACKER
OUT=$DIR/binpacker/
mkdir -p $OUT
cd $OUT
BinPacker -s fq -p pair -l $R1 -r $R2 -o $DIR/binpacker
cd ..
ln -s $OUT/BinPacker_Out_Dir/BinPacker.fa $FINAL_DIR/binpacker.fasta
###################################################################################

###################################################################################
## SHANNON
mkdir $DIR/shannon
python shannon.py -p $THREADS -o $DIR/shannon --left $R1 --right $R2
ln -s $DIR/shannon/shannon.fasta $FINAL_DIR/shannon.fasta
###################################################################################

###################################################################################
## BRIDGER
OUT=$DIR/bridger/
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT # start second time, because first run crashs at some point and then finishs
ln -s $OUT/Bridger.fasta $FINAL_DIR/bridger.fasta
###################################################################################

###################################################################################
## MIRA
mkdir $DIR/mira
touch $DIR/mira/manifest.config
echo "project = $PROJECT" > $DIR/mira/manifest.config
echo "job = est,denovo,accurate" >> $DIR/mira/manifest.config
echo "parameters = -DI:trt=/mnt/dessertlocal/projects/transcriptome_assembly/tmp -NW:cnfs=no -NW:cmrnl=no -SK:mmhr=1" >> $DIR/mira/manifest.config
# if @type == :SE
# config << "# defining illumina single end reads"
# config << "readgroup = se"
# config << "data = #{@reads}"
# config << "technology = #{@technology}"
# end
echo "# defining illumina paired end reads" >> $DIR/mira/manifest.config
echo "readgroup = pe" >> $DIR/mira/manifest.config
echo "data = $R1 $R2" >> $DIR/mira/manifest.config
echo "technology = solexa" >> $DIR/mira/manifest.config
echo "template_size = 50 1000 autorefine" >> $DIR/mira/manifest.config
echo "segment_placement = ---> <---" >> $DIR/mira/manifest.config

mira -c $DIR/mira -t $THREADS $DIR/mira/manifest.config
###################################################################################



Assembly of Homo sapiens + EBOV 7h (download script)
#######################################
## EBOV HSA 7h
## paired-end, not strand specific, 100bp
#######################################

##PARAMETERS SHARED
PROJECT=ebov_hsa_7h
R1=ebov_hsa_7h_1.fastq
R2=ebov_hsa_7h_2.fastq
THREADS=48
MEM=400

DIR=~/$PROJECT
FINAL_DIR=~/$PROJECT/final

mkdir -p $DIR
mkdir -p $FINAL_DIR

kmer1=25
kmer2=29
kmer3=33
kmer4=37
kmer5=41


###################################################################################
## TRINITY
mkdir -p $DIR/trinity
Trinity --seqType fq --max_memory $MEM"G" --left $R1 --right $R2 --CPU $THREADS --output $DIR/trinity/
ln -s $DIR/trinity/Trinity.fasta $FINAL_DIR/trinity.fasta
###################################################################################

###################################################################################
## TRANS-ABYSS
mkdir -p $DIR/transabyss

for kmer in $kmer1 $kmer2 $kmer3 $kmer4 $kmer5; do
name=k${kmer}
assemblydir=$DIR/transabyss/${name}
mkdir -p $assemblydir
finalassembly=${assemblydir}/${name}-final.fa
transabyss -k ${kmer} --pe $R1 $R2 --outdir ${assemblydir} --name ${name} --threads $THREADS
done

mergedassembly=$DIR/transabyss/merged.fa
finalassembly1=$DIR/transabyss/k"$kmer1"/k"$kmer1"-final.fa
finalassembly2=$DIR/transabyss/k"$kmer2"/k"$kmer2"-final.fa
finalassembly3=$DIR/transabyss/k"$kmer3"/k"$kmer3"-final.fa
finalassembly4=$DIR/transabyss/k"$kmer4"/k"$kmer4"-final.fa
finalassembly5=$DIR/transabyss/k"$kmer5"/k"$kmer5"-final.fa
transabyss-merge --threads $THREADS --mink 25 --maxk 41 --prefixes k${kmer1}. k${kmer2}. k${kmer3}. k${kmer4}. k${kmer5}. --out ${mergedassembly} ${finalassembly1} ${finalassembly2} ${finalassembly3} ${finalassembly4} ${finalassembly5}
ln -s $DIR/transabyss/merged.fa $FINAL_DIR/trans-abyss.fasta
###################################################################################

###################################################################################
## SOAP-TRANS
mkdir -p $DIR/soap-trans
CONFIG=$DIR/soap-trans/test.config
touch $CONFIG
echo '#maximal read length' > $CONFIG
echo "max_rd_len=100" >> $CONFIG
echo "[LIB]" >> $CONFIG
echo "#maximal read length in this lib" >> $CONFIG
echo "rd_len_cutoff=100" >> $CONFIG
#echo "#average insert size" >> $CONFIG
#echo "avg_ins=200" >> $CONFIG
echo "#if sequence needs to be reversed" >> $CONFIG
echo "reverse_seq=0" >> $CONFIG
echo "#in which part(s) the reads are used" >> $CONFIG
echo "asm_flags=3" >> $CONFIG
#echo "q=#{@reads}" >> $CONFIG #if @type == :SE
echo "#fastq paired end files" >> $CONFIG
echo "q1="$R1 >> $CONFIG
echo "q2="$R2 >> $CONFIG

kmer_min=$kmer1
kmer_max=$kmer5

# run soap with defualt settings
mkdir -p $DIR/soap-trans/default
SOAPdenovo-Trans-127mer all -s $DIR/soap-trans/test.config -o $DIR/soap-trans/default/k23 -p $THREADS
ln -s $DIR/soap-trans/default/k23.scafSeq $FINAL_DIR/soap-trans-default.fasta
###################################################################################

###################################################################################
## RNA-SPADES
mkdir -p $DIR/rnaspades/
rnaspades.py -1 $R1 -2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/rnaspades/
ln -s $DIR/rnaspades/transcripts.fasta $FINAL_DIR/rna-spades.fasta
###################################################################################

###################################################################################
## SPADES
mkdir -p $DIR/spades
spades.py --sc -1 $R1 -2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/spades/
ln -s $DIR/spades/scaffolds.fasta $FINAL_DIR/spades.fasta
###################################################################################

###################################################################################
## OASES
k_min=$kmer1
k_max=$((kmer5 + 1))
k_step=4
type=-shortPaired

mkdir -p $DIR/oases
oases_pipeline.py -m $k_min -M $k_max -s $k_step -o $DIR/oases/ -d " $type -fastq $R1 $R2 -separate " -c
ln -s $DIR/oases/Merged/transcripts.fa $FINAL_DIR/oases.fasta
###################################################################################

###################################################################################
## IDBA-TRANS
###this tools assume the paired-end reads are in order (->, <-). If your data is in order (<-, ->), please convert it by yourself.
mkdir $DIR/idba-tran
fq2fa --merge --filter $R1 $R2 $DIR/idba-tran/r12.fasta
idba_tran -r $DIR/idba-tran/r12.fasta -o $DIR/idba-tran/ --mink 25 --maxk 41 --step 4
ln -s $DIR/idba-tran/contig.fa $FINAL_DIR/idba-tran.fasta
###################################################################################

###################################################################################
## BINPACKER
OUT=$DIR/binpacker/
mkdir -p $OUT
cd $OUT
BinPacker -s fq -p pair -l $R1 -r $R2 -o $DIR/binpacker
cd ..
ln -s $OUT/BinPacker_Out_Dir/BinPacker.fa $FINAL_DIR/binpacker.fasta
###################################################################################

###################################################################################
## SHANNON
mkdir $DIR/shannon
python shannon.py -p $THREADS -o $DIR/shannon --left $R1 --right $R2
ln -s $DIR/shannon/shannon.fasta $FINAL_DIR/shannon.fasta
###################################################################################

###################################################################################
## BRIDGER
OUT=$DIR/bridger/
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT # start second time, because first run crashs at some point and then finishs
ln -s $OUT/Bridger.fasta $FINAL_DIR/bridger.fasta
###################################################################################

###################################################################################
## MIRA
mkdir $DIR/mira
touch $DIR/mira/manifest.config
echo "project = $PROJECT" > $DIR/mira/manifest.config
echo "job = est,denovo,accurate" >> $DIR/mira/manifest.config
echo "parameters = -DI:trt=/mnt/dessertlocal/projects/transcriptome_assembly/tmp -NW:cnfs=no -NW:cmrnl=no -SK:mmhr=1" >> $DIR/mira/manifest.config
# if @type == :SE
# config << "# defining illumina single end reads"
# config << "readgroup = se"
# config << "data = #{@reads}"
# config << "technology = #{@technology}"
# end
echo "# defining illumina paired end reads" >> $DIR/mira/manifest.config
echo "readgroup = pe" >> $DIR/mira/manifest.config
echo "data = $R1 $R2" >> $DIR/mira/manifest.config
echo "technology = solexa" >> $DIR/mira/manifest.config
echo "template_size = 50 1000 autorefine" >> $DIR/mira/manifest.config
echo "segment_placement = ---> <---" >> $DIR/mira/manifest.config

mira -c $DIR/mira -t $THREADS $DIR/mira/manifest.config
###################################################################################


Assembly of Homo sapiens + EBOV 23h (download script)
#######################################
## EBOV HSA 23h
## paired-end, not strand specific, 100bp
#######################################

##PARAMETERS SHARED
PROJECT=ebov_hsa_23h
R1=ebov_hsa_23h_1.fastq
R2=ebov_hsa_23h_2.fastq
THREADS=48
MEM=400
DIR=~/$PROJECT
FINAL_DIR=~/$PROJECT/final

mkdir -p $DIR
mkdir -p $FINAL_DIR

kmer1=25
kmer2=29
kmer3=33
kmer4=37
kmer5=41


###################################################################################
## TRINITY
mkdir -p $DIR/trinity
Trinity --seqType fq --max_memory $MEM"G" --left $R1 --right $R2 --CPU $THREADS --output $DIR/trinity/
ln -s $DIR/trinity/Trinity.fasta $FINAL_DIR/trinity.fasta
###################################################################################

###################################################################################
## TRANS-ABYSS
mkdir -p $DIR/transabyss

for kmer in $kmer1 $kmer2 $kmer3 $kmer4 $kmer5; do
name=k${kmer}
assemblydir=$DIR/transabyss/${name}
mkdir -p $assemblydir
finalassembly=${assemblydir}/${name}-final.fa
transabyss -k ${kmer} --pe $R1 $R2 --outdir ${assemblydir} --name ${name} --threads $THREADS
done

mergedassembly=$DIR/transabyss/merged.fa
finalassembly1=$DIR/transabyss/k"$kmer1"/k"$kmer1"-final.fa
finalassembly2=$DIR/transabyss/k"$kmer2"/k"$kmer2"-final.fa
finalassembly3=$DIR/transabyss/k"$kmer3"/k"$kmer3"-final.fa
finalassembly4=$DIR/transabyss/k"$kmer4"/k"$kmer4"-final.fa
finalassembly5=$DIR/transabyss/k"$kmer5"/k"$kmer5"-final.fa
transabyss-merge --threads $THREADS --mink 25 --maxk 41 --prefixes k${kmer1}. k${kmer2}. k${kmer3}. k${kmer4}. k${kmer5}. --out ${mergedassembly} ${finalassembly1} ${finalassembly2} ${finalassembly3} ${finalassembly4} ${finalassembly5}
ln -s $DIR/transabyss/merged.fa $FINAL_DIR/trans-abyss.fasta
###################################################################################

###################################################################################
## SOAP-TRANS
mkdir -p $DIR/soap-trans
CONFIG=$DIR/soap-trans/test.config
touch $CONFIG
echo '#maximal read length' > $CONFIG
echo "max_rd_len=100" >> $CONFIG
echo "[LIB]" >> $CONFIG
echo "#maximal read length in this lib" >> $CONFIG
echo "rd_len_cutoff=100" >> $CONFIG
#echo "#average insert size" >> $CONFIG
#echo "avg_ins=200" >> $CONFIG
echo "#if sequence needs to be reversed" >> $CONFIG
echo "reverse_seq=0" >> $CONFIG
echo "#in which part(s) the reads are used" >> $CONFIG
echo "asm_flags=3" >> $CONFIG
#echo "q=#{@reads}" >> $CONFIG #if @type == :SE
echo "#fastq paired end files" >> $CONFIG
echo "q1="$R1 >> $CONFIG
echo "q2="$R2 >> $CONFIG

kmer_min=$kmer1
kmer_max=$kmer5

# run soap with defualt settings
mkdir -p $DIR/soap-trans/default
SOAPdenovo-Trans-127mer all -s $DIR/soap-trans/test.config -o $DIR/soap-trans/default/k23 -p $THREADS
ln -s $DIR/soap-trans/default/k23.scafSeq $FINAL_DIR/soap-trans-default.fasta
###################################################################################

###################################################################################
## RNA-SPADES
mkdir -p $DIR/rnaspades/
rnaspades.py -1 $R1 -2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/rnaspades/
ln -s $DIR/rnaspades/transcripts.fasta $FINAL_DIR/rna-spades.fasta
###################################################################################

###################################################################################
## SPADES
mkdir -p $DIR/spades
spades.py --sc -1 $R1 -2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/spades/
ln -s $DIR/spades/scaffolds.fasta $FINAL_DIR/spades.fasta
###################################################################################

###################################################################################
## OASES
k_min=$kmer1
k_max=$((kmer5 + 1))
k_step=4
type=-shortPaired

mkdir -p $DIR/oases
oases_pipeline.py -m $k_min -M $k_max -s $k_step -o $DIR/oases/ -d " $type -fastq $R1 $R2 -separate " -c
ln -s $DIR/oases/Merged/transcripts.fa $FINAL_DIR/oases.fasta
###################################################################################

###################################################################################
## IDBA-TRANS
###this tools assume the paired-end reads are in order (->, <-). If your data is in order (<-, ->), please convert it by yourself.
mkdir $DIR/idba-tran
fq2fa --merge --filter $R1 $R2 $DIR/idba-tran/r12.fasta
idba_tran -r $DIR/idba-tran/r12.fasta -o $DIR/idba-tran/ --mink 25 --maxk 41 --step 4
ln -s $DIR/idba-tran/contig.fa $FINAL_DIR/idba-tran.fasta
###################################################################################

###################################################################################
## BINPACKER
OUT=$DIR/binpacker/
mkdir -p $OUT
cd $OUT
BinPacker -s fq -p pair -l $R1 -r $R2 -o $DIR/binpacker
cd ..
ln -s $OUT/BinPacker_Out_Dir/BinPacker.fa $FINAL_DIR/binpacker.fasta
###################################################################################

###################################################################################
## SHANNON
mkdir $DIR/shannon
python shannon.py -p $THREADS -o $DIR/shannon --left $R1 --right $R2
ln -s $DIR/shannon/shannon.fasta $FINAL_DIR/shannon.fasta
###################################################################################

###################################################################################
## BRIDGER
OUT=$DIR/bridger/
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT # start second time, because first run crashs at some point and then finishs
ln -s $OUT/Bridger.fasta $FINAL_DIR/bridger.fasta
###################################################################################

###################################################################################
## MIRA
mkdir $DIR/mira
touch $DIR/mira/manifest.config
echo "project = $PROJECT" > $DIR/mira/manifest.config
echo "job = est,denovo,accurate" >> $DIR/mira/manifest.config
echo "parameters = -DI:trt=/mnt/dessertlocal/projects/transcriptome_assembly/tmp -NW:cnfs=no -NW:cmrnl=no -SK:mmhr=1" >> $DIR/mira/manifest.config
# if @type == :SE
# config << "# defining illumina single end reads"
# config << "readgroup = se"
# config << "data = #{@reads}"
# config << "technology = #{@technology}"
# end
echo "# defining illumina paired end reads" >> $DIR/mira/manifest.config
echo "readgroup = pe" >> $DIR/mira/manifest.config
echo "data = $R1 $R2" >> $DIR/mira/manifest.config
echo "technology = solexa" >> $DIR/mira/manifest.config
echo "template_size = 50 1000 autorefine" >> $DIR/mira/manifest.config
echo "segment_placement = ---> <---" >> $DIR/mira/manifest.config

mira -c $DIR/mira -t $THREADS $DIR/mira/manifest.config
###################################################################################


Assembly of Homo sapiens simulated (download script)
#######################################
## HSA FLUX SIMULATED, 2x100bp, not strand specific
#######################################

##PARAMETERS SHARED
PROJECT=hsa_flux
R1=hsa_flux_1.fastq
R2=hsa_flux_2.fastq
THREADS=48
MEM=400
DIR=~/$PROJECT
FINAL_DIR=~/$PROJECT/final

mkdir -p $DIR
mkdir -p $FINAL_DIR

kmer1=25
kmer2=35
kmer3=45
kmer4=55
kmer5=65


###################################################################################
## TRINITY
mkdir -p $DIR/trinity
Trinity --seqType fq --max_memory $MEM"G" --left $R1 --right $R2 --CPU $THREADS --output $DIR/trinity/
ln -s $DIR/trinity/Trinity.fasta $FINAL_DIR/trinity.fasta
###################################################################################

###################################################################################
## TRANS-ABYSS
mkdir -p $DIR/transabyss

for kmer in $kmer1 $kmer2 $kmer3 $kmer4 $kmer5; do
name=k${kmer}
assemblydir=$DIR/transabyss/${name}
mkdir -p $assemblydir
finalassembly=${assemblydir}/${name}-final.fa
transabyss -k ${kmer} --pe $R1 $R2 --outdir ${assemblydir} --name ${name} --threads $THREADS
done

mergedassembly=$DIR/transabyss/merged.fa
finalassembly1=$DIR/transabyss/k"$kmer1"/k"$kmer1"-final.fa
finalassembly2=$DIR/transabyss/k"$kmer2"/k"$kmer2"-final.fa
finalassembly3=$DIR/transabyss/k"$kmer3"/k"$kmer3"-final.fa
finalassembly4=$DIR/transabyss/k"$kmer4"/k"$kmer4"-final.fa
finalassembly5=$DIR/transabyss/k"$kmer5"/k"$kmer5"-final.fa
transabyss-merge --threads $THREADS --mink $kmer1 --maxk $kmer5 --prefixes k${kmer1}. k${kmer2}. k${kmer3}. k${kmer4}. k${kmer5}. --out ${mergedassembly} ${finalassembly1} ${finalassembly2} ${finalassembly3} ${finalassembly4} ${finalassembly5}
ln -s $DIR/transabyss/merged.fa $FINAL_DIR/trans-abyss.fasta
###################################################################################

###################################################################################
## SOAP-TRANS
mkdir -p $DIR/soap-trans
CONFIG=$DIR/soap-trans/$PROJECT.config
touch $CONFIG
echo '#maximal read length' > $CONFIG
echo "max_rd_len=100" >> $CONFIG
echo "[LIB]" >> $CONFIG
echo "#maximal read length in this lib" >> $CONFIG
echo "rd_len_cutoff=100" >> $CONFIG
#echo "#average insert size" >> $CONFIG
#echo "avg_ins=200" >> $CONFIG
echo "#if sequence needs to be reversed" >> $CONFIG
echo "reverse_seq=0" >> $CONFIG
echo "#in which part(s) the reads are used" >> $CONFIG
echo "asm_flags=3" >> $CONFIG
#echo "q=#{@reads}" >> $CONFIG #if @type == :SE
echo "#fastq paired end files" >> $CONFIG
echo "q1="$R1 >> $CONFIG
echo "q2="$R2 >> $CONFIG

kmer_min=$kmer1
kmer_max=$kmer5

# run soap with defualt settings
mkdir -p $DIR/soap-trans/default
SOAPdenovo-Trans-127mer all -s $DIR/soap-trans/$PROJECT.config -o $DIR/soap-trans/default/k23 -p $THREADS
ln -s $DIR/soap-trans/default/k23.scafSeq $FINAL_DIR/soap-trans-default.fasta
###################################################################################

###################################################################################
## RNA-SPADES
mkdir -p $DIR/rnaspades
rnaspades.py -1 $R1 -2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/rnaspades/
ln -s $DIR/rnaspades/transcripts.fasta $FINAL_DIR/rna-spades.fasta
###################################################################################

###################################################################################
## SPADES
mkdir -p $DIR/spades
spades.py --sc -1 $R1 -2 $R2 -t $THREADS -m $MEM --cov-cutoff auto -o $DIR/spades/
ln -s $DIR/spades/scaffolds.fasta $FINAL_DIR/spades.fasta
###################################################################################

###################################################################################
## OASES
k_min=$kmer1
k_max=$((kmer5 + 1))
k_step=10
type=-shortPaired

mkdir -p $DIR/oases
oases_pipeline.py -m $k_min -M $k_max -s $k_step -o $DIR/oases/ -d " $type -fastq $R1 $R2 -separate " -c
ln -s $DIR/oases/Merged/transcripts.fa $FINAL_DIR/oases.fasta
###################################################################################

###################################################################################
## IDBA-TRANS
###this tools assume the paired-end reads are in order (->, <-). If your data is in order (<-, ->), please convert it by yourself.
mkdir $DIR/idba-tran
fq2fa --merge --filter $R1 $R2 $DIR/idba-tran/r12.fasta
idba_tran -r $DIR/idba-tran/r12.fasta -o $DIR/idba-tran/ --mink $kmer1 --maxk $kmer5 --step 10 --num_threads $THREADS
ln -s $DIR/idba-tran/contig.fa $FINAL_DIR/idba-tran.fasta
###################################################################################

###################################################################################
## BRIDGER
OUT=$DIR/bridger/
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT
Bridger.pl --seqType fq --left $R1 --right $R2 --CPU $THREADS -o $OUT # start second time, because first run crashs at some point and then finishs
ln -s $OUT/Bridger.fasta $FINAL_DIR/bridger.fasta
###################################################################################

###################################################################################
## BINPACKER
OUT=$DIR/binpacker/
mkdir -p $OUT
BinPacker -s fq -p pair -l $R1 -r $R2 -o $OUT
ln -s $OUT/BinPacker_Out_Dir/BinPacker.fa $FINAL_DIR/binpacker.fasta
###################################################################################

###################################################################################
## SHANNON
mkdir $DIR/shannon
python shannon.py -p $THREADS -o $DIR/shannon --left $R1 --right $R2
ln -s $DIR/shannon/shannon.fasta $FINAL_DIR/shannon.fasta
###################################################################################

###################################################################################
## MIRA
mkdir $DIR/mira
touch $DIR/mira/manifest.config
echo "project = $PROJECT" > $DIR/mira/manifest.config
echo "job = est,denovo,accurate" >> $DIR/mira/manifest.config
echo "parameters = -DI:trt=/mnt/dessertlocal/projects/transcriptome_assembly/tmp -NW:cnfs=no -NW:cmrnl=no -SK:mmhr=1" >> $DIR/mira/manifest.config
# if @type == :SE
# config << "# defining illumina single end reads"
# config << "readgroup = se"
# config << "data = #{@reads}"
# config << "technology = #{@technology}"
# end
echo "# defining illumina paired end reads" >> $DIR/mira/manifest.config
echo "readgroup = pe" >> $DIR/mira/manifest.config
echo "data = $R1 $R2" >> $DIR/mira/manifest.config
echo "technology = solexa" >> $DIR/mira/manifest.config
echo "template_size = 50 1000 autorefine" >> $DIR/mira/manifest.config
echo "segment_placement = ---> <---" >> $DIR/mira/manifest.config

mira -c $DIR/mira -t $THREADS $DIR/mira/manifest.config
ln -s $DIR/mira/"$PROJECT"_assembly/"$PROJECT"_d_results/"$PROJECT"_out.unpadded.fasta $FINAL_DIR/mira.fasta
###################################################################################












S4: HISAT2 re-mapping rate

back to top

We mapped the quality-controlled reads back to their assembled contigs, to determine the amount of reads that were included in the assembly process by each tool.

S5: RNAQuast statistics

back to top

Overview of RNA-Quast statistical output. Please click on the corresponding link to get a short summary of the assembly comparisons or to observe the full output.




S6: TransRate

back to top

TransRate: reference-free quality assessment of de novo transcriptome assemblies (27252236)




S7: ExN50

back to top

An alternative to the standard contig Nx (e.g. N50) statistic that can be more appropriate for transcriptome assembly is the ExN50 statistic. The N50 statistic is computed as usual but limited to the top most highly expressed transcripts that represent x% of the total normalized expression data. Therefore, we first calculated transcript abundances using the alignment-free quantification tool Salmon and subsequently used the Trinity utilities to calculate ExN50 values. We report the Ex90N50 metric as one of our main metrics.