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
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.