#Aligning with hisat2
python hisat2-2.1.0/extract_splice_sites.py annotation.gtf>annotation.ss
python hisat2-2.1.0/extract_exons.py annotation.gtf>annotation.exon
hisat2-2.1.0/hisat2-build -p 30 --ss annotation.ss --exon annotation.exon genome.fasta index

for id in `cat ID.txt`
do
{
#adapter removing
trim_galore --paired -o remove-noadapter tj-${id}_1.fq tj-${id}_2.fq
#Cleaning
SolexaQA/SolexaQA dynamictrim tj-${id}_1_val_1.fq
SolexaQA/SolexaQA dynamictrim tj-${id}_2_val_2.fq
SolexaQA/SolexaQA lengthsort tj-${id}_1_val_1.fq.trimmed tj-${id}_2_val_2.fq.trimmed
hisat2-2.1.0/hisat2 -x index -1 tj-${id}_1_val_1.fq.trimmed.paired -2 tj-${id}_2_val_2.fq.trimmed.paired -S ${id}.sam -p 30
#Samtools view format conversion
samtools-1.8/samtools view -bS ${id}.sam >${id}.bam
samtools-1.8/samtools sort ${id}.bam -o ${id}.bam-sort.bam -m 10G
#Quantification and expression
stringtie/stringtie ${id}.bam-sort.bam -o outRes${id}.gtf -p 28 -G ../annotation.gtf -A gene_abund${id}.tab -B -e
cufflinks/gffread -w ${id}transcript.fa -g genome.fasta outRes${id}.gtf -F
}
done

R CMD BATCH transcriptome.R