Content Comparison

...

for r1R1 in Raw_Reads/*_R1_*001.fastq.gz; do r2R2=$(echo $r1$R1| sed 's/_R1_/_R2_/'); name=$(echo $r1$R1|sed 's/_R1_001.fastq.gz//'|sed 's/Raw_Reads\///'); echo "fastp PE $r1 $r2-i $R1 -I $R2 -baseouto Trim_Reads/$name.trim.R1.fastq.gz ILLUMINACLIP:TruSeq3-PE-2.fa:4:30:10 MINLEN:30";done > trim_commands

for R1 in Raw_Reads/*_R1_001.fastq.gz; do R2=$(echo $R1| sed 's/_R1_/_R2_/'); name=$(echo $R1|sed 's/_R1_001.fastq.gz//'|sed 's/Raw_Reads\///'); echo "fastp -i $R1 -I $R2 -o Trim_Reads/$name.trim.R1.fastq.gz -O Trim_Reads/$name.trim.R2.fastq.gz -w 4 --detect_adapter_for_pe -j 04_Trim_Logs/$name.json -h 04_Trim_Logs/$name.html &> 04_Trim_Logs/$name.log.txt";done > fastp.commands

Use the head and wc -l to see what the output is and how much there is of it respectively.

Again as we discussed in the multiqc tutorial, running this number of commands is kind of a boarder line case, there are not a lot of total reads, but there are a large number of samples and our command does request additional processors so we should not be on the head node. 

-O Trim_Reads/$name.trim.R2.fastq.gz -w 4 --detect_adapter_for_pe -j Trim_Logs/$name.json -h Trim_Logs/$name.html &> Trim_Logs/$name.log.txt";done > fastp.commands

Use the head and wc -l to see what the output is and how much there is of it respectively.

Again as we discussed in the multiqc tutorial, running this number of commands is kind of a boarder line case of being better to run in idev or as a submitted job, there are not a lot of total reads, but there are a large number of samples. Because our command does request additional processors we should not run on the head node, if we wanted to use only a single processor, the job would take so long to run on the head node that even waiting in the queue system would be faster.

Do only one of the following, but do read through both options as there are different discussions about the process in each.

Submitting as a job

cp /corral-repl/utexas/BioITeam/gva_course/GVA2022.launcher.slurm trim.slurm
nano trim.slurm

Again while in nano you will edit most of the same lines you edited in the in the breseq tutorial. Note that most of these lines have additional text to the right of the line. This commented text is present to help remind you what goes on each line, leaving it alone will not hurt anything, removing it may make it more difficult for you to remember what the purpose of the line is

The changes to lines 22 and 23 are optional but will give you an idea of what types of email you could expect from TACC if you choose to use these options. Just be sure to pay attention to these 2 lines starting with a single # symbol after editing them.

Again use ctl-o and ctl-x to save the file and exit.

sbatch trim.slurm

...

Line 17 being set to -n 17 allows 17 jobs to run at the same time, since our command uses -w 4 (4 threads) this job will use all 68 threads available. The changes to lines 22 and 23 are optional but will give you an idea of what types of email you could expect from TACC if you choose to use these options. Just be sure to pay attention to these 2 lines starting with a single # symbol after editing them.

Again use ctl-o and ctl-x to save the file and exit.

mkdir Trim_Reads Trim_Logs
sbatch trim.slurm

The job should take less than 10 minutes once it starts if everything is working correctly, the showq -u command can be used to check for the job finishing.

Running on idev

The alternative to running all the commands as a submitted job is of course to run on an idev node. In the multiqc tutorial, we were able to tell fastqc to analyze all samples just by using the "*" wildcard as the only required input to fastqc is the filename. Here our command is much more intricate which may seem like it precludes us from being able to run interactively as we never would type 272 nearly identical commands. Obviously there is a trick for this.

chmod +x <filename>
./<filename>

Before we jump to making our commands file executable and executing it, we want to change it to be slightly different. Specifically, above we used -w 4 to specify we wanted to use 4 processor for every command. While this worked great when we also were launching 17 processes at the same time as it used all 68 processes, when executing a commands file from the command line without the help of the queue system, only 1 sample at a time will launch so you likely think we need to increase to 68 processors. 

#1 change the for loop:
for R1 in Raw_Reads/*_R1_001.fastq.gz; do R2=$(echo $R1| sed 's/_R1_/_R2_/'); name=$(echo $R1|sed 's/_R1_001.fastq.gz//'|sed 's/Raw_Reads\///'); echo "fastp -i $R1 -I $R2 -o Trim_Reads/$name.trim.R1.fastq.gz -O Trim_Reads/$name.trim.R2.fastq.gz -w 68 --detect_adapter_for_pe -j Trim_Logs/$name.json -h Trim_Logs/$name.html &> Trim_Logs/$name.log.txt";done > fastp.commands

#2 use sed to do an in file replacement (something new you haven't seen before in this class)

14 minutes with -w4 idev

16 minutes with -w 68 idev

error after 54 samples -w 68 idev & ... <1 minute

error after 148 samples -w 1 idev & ... <1 minute

~2 minutes as submitted job -w4 -n 17

? minutes as submitted job with -w68 -n 1

? minutes as submitted job with -w1 -n 68

Evaluating the output

ls Trim_Reads | wc -l
grep -c "TrimmomaticPE: Completed successfully" Queue_job.o*
grep -c "100.00%" Queue_job.o*

...