1.filter_linker.pbs

#!/bin/bash
#PBS -l nodes=1:ppn=20
#PBS -l walltime=12:00:00
#PBS -l mem=60GB
#PBS -l vmem=60GB
#PBS -j oe


# ChIA-PIPE
#         Step 1: Filter the linker sequence
#         (Requires single-linker ChIA-PET data)
# 2018
# The Jackson Laboratory for Genomic Medicine


## The help message:
function usage
{
    echo -e "usage: qsub -F \"--conf ${conf} --out_dir ${out_dir}\" 
    1.filter_linker.pbs
    " 
}

## Parse arguments from the command line
while [ "$1" != "" ]; do
    case $1 in
        -c | --conf )           shift
                                conf=$1
                                ;;
        -o | --out_dir )        shift
                                out_dir=$1
                                ;;
        -h | --help )           usage
                                exit
                                ;;
        * )                     usage
                                exit 1
    esac
    shift
done


## Move to output directory
cd ${PBS_O_WORKDIR}
source ${conf}
mkdir -p ${out_dir}
cd ${out_dir}

# Add dependency dir to path
export PATH=${dep_dir}:${PATH}

## Create name of the log file
log_file=1.${run}.filter_linker.log

## Report arguments parsed
echo "
Arguments:
    main_prog=${main_prog}
    run=${run}
    data_dir=${data_dir}
    out_dir=${out_dir}
" >> ${log_file}


## Perform linker detection and generate different categories of fastq files
# Report linker detection start
echo "
`date` --- Linker detection started on: ---
    ${data_dir}/${r1_fastq},
    ${data_dir}/${r2_fastq}
" >> ${log_file}


## Linker filtering
if [ ${experiment_type} == 'ChIA-PET' ]; then
    # Linker filtering for ChIA-PET experiment

    if [ ${linker_a} == "none" ]; then
        # Default linker sequence
        # Single linker with default sequence
        ${main_prog} stag -W -T ${min_tag_len} -t ${n_thread} -O ${run} \
            ${data_dir}/${r1_fastq} ${data_dir}/${r2_fastq} \
            2>> ${log_file}
    
    elif [ ${linker_a} != "none" ]; then
        # Custom linker sequence
        if [ ${linker_b} == "none" ]; then
            # Single linker with custom sequence
            ${main_prog} stag -W -T ${min_tag_len} -t ${n_thread} -O ${run} \
                -A ${linker_a} \
                ${data_dir}/${r1_fastq} ${data_dir}/${r2_fastq} \
                2>> ${log_file}
        elif [ ${linker_b} != "none" ]; then
            # Two linkers with custom sequence
                ${main_prog} tag -W -T ${min_tag_len} -t ${n_thread} -O ${run} \
                -A ${linker_a} -B ${linker_b}\
                ${data_dir}/${r1_fastq} ${data_dir}/${r2_fastq} \
                2>> ${log_file}
        fi
    fi

    # Report linker detection completion
    echo -e "`date` --- Linker detection completed ----\n" >> ${log_file}


    ## Get the statistics
    # Report statistics start
    echo "
    `date` --- Statistics started  ----
    " >> ${log_file}


    if [ ${linker_b} == "none" ]; then
        # Statistics
        ${main_prog} stat -s -p -T ${min_tag_len} -t ${n_thread} ${run}.cpu \
            2>> ${log_file} 1> ${run}.stat
    else
        ${main_prog} stat -d -p -T ${min_tag_len} -t ${n_thread} ${run}.cpu \
            2>> ${log_file} 1> ${run}.stat    
    fi

    # Report statistics completion
    echo -e "`date` --- Statistics completed  ----\n" >> ${log_file}

elif [ ${experiment_type} == 'HiChIP' ]; then
    # Linker filtering for HiChIP data
    ${dep_dir}/python ${bin_dir}/util/scripts/filter_hichip_linker.py \
        --r1_file  ${data_dir}/${r1_fastq} \
        --r2_file  ${data_dir}/${r2_fastq} \
        --run  ${run} \
        --linker  ${linker_a} \
        --min_tag_len  ${min_tag_len}
    
    # Write linker filtering stats
    bash ${bin_dir}/util/scripts/write_hichip_linker_stats.sh -c ${conf}

elif [ ${experiment_type} == 'PLAC-seq' ]; then
    # Linker filtering for PLAC-seq data
    ${dep_dir}/python ${bin_dir}/util/scripts/filter_hichip_linker.py \
        --r1_file  ${data_dir}/${r1_fastq} \
        --r2_file  ${data_dir}/${r2_fastq} \
        --run  ${run} \
        --linker  ${linker_a} \
        --min_tag_len  ${min_tag_len}
    
    # Write linker filtering stats
    bash ${bin_dir}/util/scripts/write_hichip_linker_stats.sh -c ${conf}
    
fi

## Compress files
pigz -p ${n_thread} ${run}.singlelinker.paired.fastq 2>> ${log_file}
pigz -p ${n_thread} ${run}.none.fastq 2>> ${log_file}

pigz -p ${n_thread} ${run}.singlelinker.single.fastq 2>> ${log_file}
pigz -p ${n_thread} ${run}.conflict.fastq 2>> ${log_file}
pigz -p ${n_thread} ${run}.tied.fastq 2>> ${log_file}


## Report script completion
echo -e "$0 done \n" >> ${log_file}
echo "`date`" >> ${log_file}


# Submit next step of pipeline
if [ ${all_steps} == true ]
then
    # Set the resource parameters for the computing cluster
    # depending on the run type (miseq or hiseq)
    
    if [ ${run_type} == "miseq" ]
    then
        wall_time=10
    elif [ ${run_type} == "hiseq" ] || [ ${run_type} == "nextseq" ]
    then
        wall_time=40
    else
        wall_time=60
    fi
    
    ### 2. Mapping
    # Submit a separate job for each category of reads
    # (none, singlelinker.single, singlelinker.paired)
    job_2a=$( qsub -F \
    "--conf ${conf} --out_dir ${out_dir} --tag_name none" \
    -l nodes=1:ppn=${n_thread},mem=${mem}gb,vmem=${mem}gb,walltime=${wall_time}:00:00 \
    -j oe -o ${out_dir}/2.${run}.map.none.o \
    ${bin_dir}/2.map.pbs )


    job_2b=$( qsub -F \
    "--conf ${conf} --out_dir ${out_dir} --tag_name singlelinker.single" \
    -l nodes=1:ppn=${n_thread},mem=${mem}gb,vmem=${mem}gb,walltime=${wall_time}:00:00 \
    -j oe -o ${out_dir}/2.${run}.map.singlelinker.single.o \
    ${bin_dir}/2.map.pbs )


    job_2c=$( qsub -F \
    "--conf ${conf} --out_dir ${out_dir} --tag_name singlelinker.paired" \
    -l nodes=1:ppn=${n_thread},mem=${mem}gb,vmem=${mem}gb,walltime=${wall_time}:00:00 \
    -j oe -o ${out_dir}/2.${run}.map.singlelinker.paired.o \
    ${bin_dir}/2.map.pbs )

fi