Skip to content

Commit

Permalink
Switched from bwa aln to mem.
Browse files Browse the repository at this point in the history
Not fully switched for single end reads.
  • Loading branch information
bergeycm committed Feb 19, 2016
1 parent d864e32 commit e7e8382
Show file tree
Hide file tree
Showing 4 changed files with 13 additions and 50 deletions.
37 changes: 7 additions & 30 deletions full_analysis.mk
Expand Up @@ -27,11 +27,9 @@ else ifeq ($(READ_TYPE),PE)
endif
# --- alignment_steps
ifeq ($(READ_TYPE),SE)
align : results/${IND_ID}.readSE.bwa.${GENOME_NAME}.sai
sampe_or_samse : results/${IND_ID}.SE.bwa.${GENOME_NAME}.sam
align : results/${IND_ID}.SE.bwa.${GENOME_NAME}.sam
else ifeq ($(READ_TYPE),PE)
align : results/${IND_ID}.read1.bwa.${GENOME_NAME}.sai
sampe_or_samse : results/${IND_ID}.PE.bwa.${GENOME_NAME}.sam
align : results/${IND_ID}.PE.bwa.${GENOME_NAME}.sam
endif
sam2bam : results/${IND_ID_W_PE_SE}.bwa.${GENOME_NAME}.sam.bam
sort_and_index_bam : results/${IND_ID_W_PE_SE}.bwa.${GENOME_NAME}.sam.bam.sorted.bam.bai
Expand All @@ -56,7 +54,7 @@ compress_and_upload : results/${IND_ID_W_PE_SE}.bwa.${GENOME_NAME}.tar.gz
# --- Individual steps
preliminary_steps : index_genome make_dict
pre_aln_analysis_steps : fastqc
alignment_steps : align sampe_or_samse sam2bam sort_and_index_bam get_alignment_stats
alignment_steps : align sam2bam sort_and_index_bam get_alignment_stats
post_alignment_filtering_steps : fix_mate_pairs filter_unmapped add_read_groups filter_bad_qual
snp_calling_steps : local_realign_targets local_realign call_snps filter_snps get_snp_stats call_consensus
archive_steps : compress_and_upload
Expand Down Expand Up @@ -148,38 +146,17 @@ reports/${IND_ID}.readSE.stats.zip : ${READ_SE} ${FASTQC}/* #scripts/run_fastqc.
# --- Align reads to genome with BWA
# -------------------------------------------------------------------------------------- #

# Alignment output (*.sai) depends on bwa, the reads FASTAs, the genome (index), and align.sh
# Using the first read as a stand in for the both
results/${IND_ID}.read1.bwa.${GENOME_NAME}.sai : ${BWA}/* ${READ1} ${READ2} ${GENOME_FA}i #scripts/align.sh
# Alignment output (*.sam) depends on bwa, the reads FASTAs, the genome (index), and align.sh
results/${IND_ID}.PE.bwa.${GENOME_NAME}.sam : ${BWA}/* ${READ1} ${READ2} ${GENOME_FA}i #scripts/align.sh
@echo "# === Aligning reads to genome ================================================ #";
./scripts/align.sh ${GENOME_FA} ${GENOME_NAME};

# Read 2 depends on read 1
ifeq ($(READ_TYPE),PE)
results/${IND_ID}.read2.bwa.${GENOME_NAME}.sai : results/${IND_ID}.read1.bwa.${GENOME_NAME}.sai
endif

# Align SE reads
# Alignment output (*.sai) depends on bwa, the reads FASTAs, the genome (index), and alignSE.sh
results/${IND_ID}.readSE.bwa.${GENOME_NAME}.sai : ${BWA}/* ${READ_SE} ${GENOME_FA}i #scripts/align.sh
# Alignment output (*.sam) depends on bwa, the reads FASTAs, the genome (index), and alignSE.sh
results/${IND_ID}.SE.bwa.${GENOME_NAME}.sam : ${BWA}/* ${READ_SE} ${GENOME_FA}i #scripts/align.sh
@echo "# === Aligning SE reads to genome ============================================= #";
./scripts/alignSE.sh ${GENOME_FA} ${GENOME_NAME};

# -------------------------------------------------------------------------------------- #
# --- Run sampe and samse to generate SAM files
# -------------------------------------------------------------------------------------- #

# sampe output (*.sam) depends on *.sai files and sampe.sh
# Using the first read as a stand in for the both
results/${IND_ID}.PE.bwa.${GENOME_NAME}.sam : results/${IND_ID}.read1.bwa.${GENOME_NAME}.sai #scripts/sampe.sh
@echo "# === Making SAM file from PE reads =========================================== #";
./scripts/sampe.sh ${GENOME_FA} ${GENOME_NAME};

# samse output (*.sam) depends on *.sai file and sampe.sh
results/${IND_ID}.SE.bwa.${GENOME_NAME}.sam : results/${IND_ID}.readSE.bwa.${GENOME_NAME}.sai #scripts/sampe.sh
@echo "# === Making SAM file from SE reads =========================================== #";
./scripts/samse.sh ${GENOME_FA} ${GENOME_NAME};

# -------------------------------------------------------------------------------------- #
# --- Convert SAM file to BAM file
# -------------------------------------------------------------------------------------- #
Expand Down
2 changes: 1 addition & 1 deletion indiv_analysis
Expand Up @@ -4,4 +4,4 @@
# This compares a single individual
# Input data must be specified by editing config.mk

make -s -f full_analysis.mk indiv
make -s -f full_analysis.mk indiv
23 changes: 4 additions & 19 deletions scripts/align.sh
@@ -1,7 +1,7 @@
#!/bin/sh

# ------------------------------------------------------------------------------
# --- Align reads to genome with BWA
# --- Align reads to genome with BWA mem
# ------------------------------------------------------------------------------

# Check that genome FASTA, and genome code were passed as parameters
Expand All @@ -15,24 +15,9 @@ fi
GENOME_PATH=$(echo $1 | sed 's/.[^.]*$//g')
GENOME_CODE=$2

echo "CMD: $BWA/bwa aln \
$BWA/bwa mem \
$BWA_ALN_PARAM \
$GENOME_PATH \
$READ1 > results/${IND_ID}.read1.bwa.${GENOME_CODE}.sai";
$READ1 $READ2 > results/${IND_ID}.PE.bwa.${GENOME_CODE}.sam

$BWA/bwa aln \
$BWA_ALN_PARAM \
$GENOME_PATH \
$READ1 > results/${IND_ID}.read1.bwa.${GENOME_CODE}.sai

echo "CMD: $BWA/bwa aln \
$BWA_ALN_PARAM \
$GENOME_PATH \
$READ2 > results/${IND_ID}.read2.bwa.${GENOME_CODE}.sai";

$BWA/bwa aln \
$BWA_ALN_PARAM \
$GENOME_PATH \
$READ2 > results/${IND_ID}.read2.bwa.${GENOME_CODE}.sai

exit;
exit;
1 change: 1 addition & 0 deletions test.sh
@@ -0,0 +1 @@
make --just-print -s -f full_analysis.mk indiv IND_ID=ssese1

0 comments on commit e7e8382

Please sign in to comment.