Skip to content

Commit

Permalink
fine tune the wrapper scripts
Browse files Browse the repository at this point in the history
  • Loading branch information
lh3 committed Nov 18, 2014
1 parent 80e4ecf commit 61ab329
Show file tree
Hide file tree
Showing 4 changed files with 17 additions and 15 deletions.
6 changes: 3 additions & 3 deletions NEWS.md
Original file line number Diff line number Diff line change
Expand Up @@ -9,9 +9,9 @@ mapping, we start to provide a wrapper script and precompiled binaries since
this release. Please check README-alt.md for details.

Another major addition to BWA-MEM is HLA typing, which made possible with the
new ALT mapping strategy. Necessary data and programs are also included in
the binary release. The wrapper script also performs HLA typing when HLA genes
are also included in the reference genome as additional ALT contigs.
new ALT mapping strategy. Necessary data and programs are included in the
binary release. The wrapper script also performs HLA typing when HLA genes are
also included in the reference genome as additional ALT contigs.

Other notable changes to BWA-MEM:

Expand Down
7 changes: 4 additions & 3 deletions bwa.1
Original file line number Diff line number Diff line change
@@ -1,4 +1,4 @@
.TH bwa 1 "19 May 2014" "bwa-0.7.9-r783" "Bioinformatics tools"
.TH bwa 1 "18 November 2014" "bwa-0.7.11-r999" "Bioinformatics tools"
.SH NAME
.PP
bwa - Burrows-Wheeler Alignment Tool
Expand Down Expand Up @@ -233,8 +233,9 @@ more aggressive read pair. [17]
.B INPUT/OUTPUT OPTIONS:
.TP
.B -p
Assume the first input query file is interleaved paired-end FASTA/Q. See the
command description for details.
Smart pairing. If two adjacent reads have the same name, they are considered
to form a read pair. This way, paired-end and single-end reads can be mixed
in a single FASTA/Q stream.
.TP
.BI -R \ STR
Complete read group header line. '\\t' can be used in
Expand Down
17 changes: 9 additions & 8 deletions extras/run-bwamem
Original file line number Diff line number Diff line change
Expand Up @@ -4,8 +4,8 @@ use strict;
use warnings;
use Getopt::Std;

my %opts = (t=>1, n=>64);
getopts("Sadsko:R:x:t:", \%opts);
my %opts = (t=>1);
getopts("SadskHo:R:x:t:", \%opts);

die('
Usage: run-bwamem [options] <idxbase> <file1> [file2]
Expand All @@ -22,6 +22,7 @@ Options: -o STR prefix for output files [inferred from
-d mark duplicate (via samblaster)
-S for SAM/BAM input, don\'t shuffle
-s sort the output alignment (requring more RAM)
-H skip HLA typing
-k keep temporary files generated by typeHLA
Examples:
Expand All @@ -31,15 +32,15 @@ Examples:
run-bwamem -o prefix -t8 -R"@RG\tID:foo\tSM:bar" hs38d6.fa read1.fq.gz read2.fq.gz
* Remap coordinate-sorted BAM, trim Illumina PE adapters and sort the output. The BAM
may contain single-end or paired-end reads, or a mixture of the two types. The read
groups are not transferred to the output BAM.
may contain single-end or paired-end reads, or a mixture of the two types. Read groups
are not transferred to the output BAM, though.
run-bwamem -sao prefix hs38d6.fa old-srt.bam
* Remap name-grouped BAM and mark duplicates. Note that in this case, all reads from
a single library should be aligned at the same time. Paired-end only.
run-bwamem -Sdo prefix hs38d6.fa old-unsrt.bam
run-bwamem -Sdo prefix hs38d6.fa old-unsrt.bam
Output files:
Expand Down Expand Up @@ -119,7 +120,7 @@ if ($is_sam || $is_bam) {

my $bwa_opts = "-p " . ($opts{t} > 1? "-t$opts{t} " : "") . (defined($opts{x})? "-x $opts{x} " : "") . (defined($opts{R})? "-R'$opts{R}' " : "");

$cmd .= " | $root/trimadap \\\n" if defined($opts{a});
$cmd .= " | $root/trimadap 2> $prefix.log.trim \\\n" if defined($opts{a});
$cmd .= " | $root/bwa mem $bwa_opts$ARGV[0] - 2> $prefix.log.bwamem \\\n";
$cmd .= " | $root/samblaster 2> $prefix.log.dedup \\\n" if defined($opts{d});

Expand All @@ -138,8 +139,8 @@ if (-f "$ARGV[0].alt") {
my $t_sort = $opts{t} < 4? $opts{t} : 4;
$cmd .= defined($opts{s})? " | $root/samtools sort -@ $t_sort -m1G - $prefix.aln;\n" : " | $root/samtools view -1 - > $prefix.aln.bam;\n";

if ($has_hla && (!defined($opts{x}) || $opts{x} eq 'intractg')) {
$cmd .= "$root/run-HLA ". ($opts{x} eq 'intractg'? "-A " : "") . "$prefix.hla > $prefix.hla.top 2> $prefix.log.hla;\n";
if ($has_hla && !defined($opts{H}) && (!defined($opts{x}) || $opts{x} eq 'intractg')) {
$cmd .= "$root/run-HLA ". (defined($opts{x}) && $opts{x} eq 'intractg'? "-A " : "") . "$prefix.hla > $prefix.hla.top 2> $prefix.log.hla;\n";
$cmd .= "cat $prefix.hla.HLA*.gt | grep ^GT | cut -f2- > $prefix.hla.all;\n";
$cmd .= "rm -f $prefix.hla.HLA*;\n" unless defined($opts{k});
}
Expand Down
2 changes: 1 addition & 1 deletion extras/typeHLA.sh
Original file line number Diff line number Diff line change
Expand Up @@ -20,7 +20,7 @@ if [ $is_ctg -eq 0 ]; then
echo "** De novo assembling..." >&2
len=`$root/seqtk comp $pre.fq | awk '{++x;y+=$2}END{printf("%.0f\n", y/x)}'`
$root/fermi2.pl unitig -f $root/fermi2 -r $root/ropebwt2 -t2 -l$len -p $pre.tmp $pre.fq > $pre.tmp.mak
make -f $pre.tmp.mak
make -f $pre.tmp.mak >&2
cp $pre.tmp.mag.gz $pre.mag.gz
else
rm -f $pre.tmp.mag.gz
Expand Down

0 comments on commit 61ab329

Please sign in to comment.