Module 2.2 fastp

Starting on sphinx from this directory: /lustre/isaac24/proj/UTK0505/yourusesrid with conda using the snps environment.

Fastp

Get set up. make a new analysis folder for this and symbolically link the input files

cd analyses
mkdir 02_fastp
cd 02_fastp
ln -s ../../raw_data/*gz .
mkdir fastp_trimmed

We never took a peak at any of our fastq.gz files to make sure they are correctly formatted. Lets try that now (pick a sample name from your list). zcat is a command-line utility in Unix and Linux used to view the uncompressed contents of a compressed file (such as .gz files) without having to manually extract or uncompress it to disk first. Very useful.

zcat 21_FP_OH_FSB_F_10_R1.fastq.gz | head -n 16

Lets run one sample first. We can generate two types of output - html and json. If we want to use MultiQC later, we need to name the json output files matching the pattern *.fastp.json

fastp   --in1 21_FP_OH_FSB_F_10_R1.fastq.gz \
        --out1 21_FP_OH_FSB_F_10_R1.trimmed.fastq.gz \
        --in2 21_FP_OH_FSB_F_10_R2.fastq.gz \
        --out2 21_FP_OH_FSB_F_10_R2.trimmed.fastq.gz \
        --detect_adapter_for_pe \
        --thread 2 \
        --json 21_FP_OH_FSB_F_10.fastp.json \
        --html 21_FP_OH_FSB_F_10.fastp.html 

It autodetected our adapter which is convenient. Secure copy the html to your laptop to check out the results.

Now lets run all the samples. Create a run_fastp.sh file. this took about a minute per file for me.

Here is a script that I wrote. I start simple:

## Run all files through fastp
for R1 in ../../raw_data/*_R1.fastq.gz
do
    basename=$(echo "$R1" | sed 's/_R1.fastq.gz//')
    echo $R1
    echo $basename
done

Now start add the fastp part and check out those variable names

## Run all files through fastp
for R1 in ../../raw_data/*_R1.fastq.gz
do
    basename=$(echo "$R1" | xargs basename | sed 's/_R1.fastq.gz//')
    echo $R1
    echo $basename


    echo "fastp   --in1 ${basename}_R1.fastq.gz \
        --out1 ${basename}_R1.trimmed.fastq.gz \
        --in2 ${basename}_R2.fastq.gz \
        --out2 ${basename}_R2.trimmed.fastq.gz \
        --detect_adapter_for_pe \
        --thread 2 \
        --json ${basename}.fastp.json \
        --html ${basename}.fastp.html
    "
done

Now remove the echo and see if it works

## Run all files through fastp
for R1 in ../../raw_data/*_R1.fastq.gz
do
    basename=$(echo "$R1" | sed 's/_R1.fastq.gz//')
    echo $R1
    echo $basename


    fastp   --in1 ${basename}_R1.fastq.gz \
        --out1 ${basename}_R1.trimmed.fastq.gz \
        --in2 ${basename}_R2.fastq.gz \
        --out2 ${basename}_R2.trimmed.fastq.gz \
        --detect_adapter_for_pe \
        --thread 2 \
        --json ${basename}.fastp.json \
        --html ${basename}.fastp.html
    
done

But if we let AI loose we get a lot more complexity...

#!/bin/bash
# run_fastp.sh
# Run fastp on every pair of paired-end FASTQ files in this folder.
# Usage:  bash run_fastp.sh

# Stop the script immediately if any command fails.
set -e

# Make a folder to hold the fastp reports (the -p means "don't complain if it already exists").
mkdir -p fastp_reports
mkdir -p trimmed

# Loop over every R1 file in the current directory.
# The * is a wildcard: it matches ANY text, so this matches all 10 R1 files.
for R1 in *_R1.fastq.gz
do
    # Build the matching R2 filename by swapping "_R1" for "_R2" in the name.
    # The ${VAR/old/new} syntax replaces the first match of "old" with "new".
    R2="${R1/_R1/_R2}"

    # Build a short sample name by chopping "_R1.fastq.gz" off the end.
    # The ${VAR%pattern} syntax removes "pattern" from the END of the value.
    # Example: 21_FP_OH_FSB_F_10_R1.fastq.gz  ->  21_FP_OH_FSB_F_10
    SAMPLE="${R1%_R1.fastq.gz}"

    # Print a message so we can watch the progress in the terminal.
    echo "Now processing sample: $SAMPLE"

    # Safety check: make sure the R2 file actually exists before running fastp.
    if [ ! -f "$R2" ]; then
        echo "  WARNING: no R2 file found for $R1 -- skipping."
        continue    # "continue" jumps to the next sample in the loop
    fi

    # This is the fastp command, with the filenames replaced by variables.
    fastp \
        --in1 "$R1" \
        --out1 trimmed/"$SAMPLE"_R1.trimmed.fastq.gz \
        --in2 "$R2" \
        --out2 trimmed/"$SAMPLE"_R2.trimmed.fastq.gz \
        --detect_adapter_for_pe \
        --thread 2 \
        --json fastp_reports/"$SAMPLE".fastp.json \
        --html fastp_reports/"$SAMPLE".fastp.html


done

echo "All done! Reports are in the fastp_reports/ folder."

What are the positives and negatives of using AI? How could you make sure you understand what it is doing and if it is doing what you expect?

Screen

What if you want this to run and you want to walk away and return later? We can use screen for that. GNU Screen is a terminal multiplexer that lets you run a bash shell or other commands and keep them active in the background even if you disconnect.

Basic Screen Commands

• Start a new session: Type screen.
• Start a named session: Type screen -S myname
• List active sessions: Type screen -ls
• Reattach to a session: Type screen -r or screen -r myname
• Detach from inside a session: Press Ctrl + A then D.

Example for our bash script above.

screen -S fastp
bash run_fastq.sh
Ctrl + A D

Now you can exit. Come back later to see if the job finished or had any errors.

screen -r fastp

For Homework:

1. Run multiqc on the output json files from fastp and see what it says.
2. Rerun fastqc and multiqc on the trimmed files. Did fastp fix the reports? Doing this is best practice to make sure fastq fixed anything that needed fixing.