This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Hisat2 file error

Hi I am trying to run the following script but I am getting this error: gzip: : /camp/home/fernanm/working/MARTA-RNAseq-analysis/LP-files/JIM460A4_S45_L008_R1_trimmed.fastq.gz_R2_trimmed.fastq.gz: No such file or directory. To me it looks like its doubling up the end of the file so It cannot find it but I am not sure how to fix it. Thanks in advance!

###Load module
ml HISAT2/2.1.0-foss-2016b

###Store sam files in this directory
OUTPUT=$HOME/working/MARTA-RNAseq-analysis/sam
mkdir -p $OUTPUT
cd $OUTPUT

###reference genome
REF=$HOME/working/MARTA-RNAseq-analysis/reference-genome/hisat2_index/grcm38

###running HISAT2
echo "Running HISAT2"

SAMPLES=$HOME/working/MARTA-RNAseq-analysis/LP-files/*_trimmed.fastq.gz

for SAMPLE in $(ls $SAMPLES)
do
        R1=${SAMPLE}_R1_trimmed.fastq.gz
        R2=${SAMPLE}_R2_trimmed.fastq.gz
        SAM=${SAMPLE}.sam

        hisat2 -p 24 -t  -x $REF -1 $R1 -2 $R2 -S $SAM
done
rna-seq

1 answer

Yes, it doubles up because your *_trimmed.fastq.gz does call all files so both R1 and R2 per basename so the loop runs twice for each pair of files. It would be better do something like this below.

First make a file that contains all basenames, then loop through these names. Lastly I added a pipe | to directly write the hisat2 output to a bam rather than sam file. It is compressed and therefore smaller. There is these days no reason to keep sam files imho.

ls $HOME/working/MARTA-RNAseq-analysis/LP-files/*_R1_trimmed.fastq.gz \
| awk -F "_R1_trimmed" '{print $1}' > samples.txt

while read SAMPLE
  do
  R1=${SAMPLE}_R1_trimmed.fastq.gz
  R2=${SAMPLE}_R2_trimmed.fastq.gz
  BAM=${SAMPLE}.bam
  hisat2 -p 24 -t  -x $REF -1 $R1 -2 $R2 | samtools view -o $BAM
  done < samples.txt

Log in to answer this question.