I am trying to do adapter trimming, alignment and sorting for a range of large scale paired end fastq files. The code I am using is given below:
#!/bin/bash
# Do Trimming, Alignment according to Config file
#SBATCH --job-name=T_Align_single
#SBATCH --nodes=1
#SBATCH --ntasks=1
#SBATCH --cpus-per-task=32
#SBATCH --partition=batch
#SBATCH --account=markay
#SBATCH --time=72:00:00
#SBATCH --mem=16G
#SBATCH --mail-type=TIME_LIMIT_80,BEGIN,END
#SBATCH --mail-user=<email-redacted>
module load bowtie2/2.4.1
module load samtools/1.9
mkdir -p Mapping/
config=${INPUT_FILE}
file1=$(cat $config | awk '$1 ~ /fastq/ {print $1}')
file2=$(cat $config | awk '$1 ~ /fastq/ {print $2}')
file3=$(cat $config | awk '$1 ~ /index=/ {print $1}'| sed 's/ref_index=//')
#adapter trimming
trimfile1=Mapping/"${file1%.fastq}".trim.fastq
trimfile2=Mapping/"${file2%.fastq}".trim.fastq
java -jar /labs/markay/Aranyak/Adriana_intro_CnT/Cnt_Test/Test-analysis/fastq_processing/Trimmomatic-0.39/trimmomatic-0.39.jar PE -phred33 $file1 $file2 \
$trimfile1 Mapping/"${file1%.fastq}".unpaired.fastq \
$trimfile2 Mapping/"${file2%.fastq}".unpaired.fastq \
LEADING:3 TRAILING:3 SLIDINGWINDOW:4:15 \
MINLEN:20 ILLUMINACLIP: /labs/markay/Aranyak/Adriana_intro_CnT/adapter.fa:2:30:10 -threads ${SLURM_CPUS_PER_TASK}
# mapping using bowtie2
sortbam=Mapping/"${file1%_R1_001.fastq}".sort.bam
sample_ID=$(cat $config | awk '$1 ~ /fastq/ {print $4}')
bowtie2 -p ${SLURM_CPUS_PER_TASK} --local --very-sensitive-local --no-unal \
--no-mixed --no-discordant --phred33 -I 10 -X 700 \
-x $file3 \
-1 $trimfile1 -2 $trimfile2 --rg-id $sample_ID | samtools view -Sb - | samtools sort -@${SLURM_CPUS_PER_TASK} -o $sortbam
#Index the sorting by coordinate bam file
samtools index ${sortbam}
My original fastq files before adapter trimming looked like this
R3-LiD-M-LK_S16_L002_R2_001.fastq.gz
R3-LiD-M-LK_S16_L002_R1_001.fastq.gz
However after adapter trimming I removed the lane names and my file looks like this
R3-LiD-M-LK_S16_R2_001.fastq.gz
R3-LiD-M-LK_S16_R1_001.fastq.gz
The code that worked for the original file is not generating a proper bam file at the end though I am getting a list of trimmed fastq files. The error message I am getting is:
[E::hts_open_format] Failed to open file R1-LiD-H-DJ_S3.sort.rmdup.bam
samtools index: failed to open "R1-LiD-H-DJ_S3.sort.rmdup.bam": No such file or directory
[E::hts_open_format] Failed to open file R1-LiD-H-DJ_S3.sort.bam
The file 'R1-LiD-H-DJ_S3.sort.bam' does not exist
[E::hts_open_format] Failed to open file R1-LiD-H-DJ_S3.sort.rmdup.bam
The file 'R1-LiD-H-DJ_S3.sort.rmdup.bam' does not exist
R1-LiD-H-DJ_S3 DONE.
How should I change the code to generate proper alignment and bam files?
set -x
after the#SBATCH
lines and look at the job's STDERR file to see the exact commands being executed. $\endgroup$R1-LiD-H-DJ_S3.sort.bam
and a file namedR1-LiD-H-DJ_S3.sort.rmdup.bam
and neither of them exist. You also don't mentionrmdub
anywhere in your code, is this really the code you are running? $\endgroup$