initial import runRAILS_lowcontiguityseqs.sh

Author: warrenlr

RAILS_v1.2/RAILS

@@ -25,7 +25,7 @@ use Getopt::Std;
 use Net::SMTP;
 use vars qw($opt_f $opt_s $opt_d $opt_i $opt_e $opt_l $opt_a $opt_v $opt_b $opt_t $opt_p $opt_q);
 getopts('f:s:d:e:l:a:v:b:t:p:i:q:');
-my ($base_name,$frag_dist,$seqid,$insert_stdev,$min_links,$max_link_ratio,$verbose)=("",1000,0.9,1.0,1,0.0,0);
+my ($base_name,$frag_dist,$seqid,$insert_stdev,$min_links,$max_link_ratio,$verbose)=("",250,0.9,1.0,1,0.0,0);
 
 my $version = "[v1.2]";
 my $dev = "rwarren\@bcgsc.ca";

RAILS_v1.2/cobbler.pl

@@ -25,7 +25,7 @@
 use Net::SMTP;
 use vars qw($opt_f $opt_s $opt_d $opt_i $opt_v $opt_b $opt_t $opt_q);
 getopts('f:s:d:v:b:t:i:q:');
-my ($base_name,$frag_dist,$seqid,$verbose)=("",1000,0.9,0);
+my ($base_name,$frag_dist,$seqid,$verbose)=("",250,0.9,0);
 
 my $version = "[v0.3]";
 my $dev = "rwarren\@bcgsc.ca";

RAILS_v1.2/readme.md

@@ -6,6 +6,7 @@
 -------------
 
 RAILS: Radial Assembly Improvement by Long Sequence Scaffolding
+
 Cobbler: Gap-filling with long sequences
 
 
@@ -89,7 +90,7 @@ Software. doi: 10.21105/joss.00116
 -------------
 
 <pre>
-./runRAILS.sh
+./runRAILS.sh or runRAILS_lowcontiguityseqs.sh
 Usage: runRAILS.sh <FASTA assembly .fa> <FASTA long sequences .fa> <anchoring sequence length eg. 250> <min sequence identity 0.95>
 
 this pipeline will:
@@ -107,7 +108,7 @@ Usage: ./cobbler.pl [v0.3]
 -f  Assembled Sequences to further scaffold (Multi-FASTA format NO LINE BREAKS, required)
 -q  Long Sequences queried (Multi-FASTA format NO LINE BREAKS, required)
 -s  BAM file (use v0.2 for reading SAM files)
--d  Anchoring bases on contig edges (ie. minimum required alignment size on contigs, default -d 1000, optional)
+-d  Anchoring bases on contig edges (ie. minimum required alignment size on contigs, default -d 250, optional)
 -i  Minimum sequence identity, default -i 0.9, optional
 -t  LIST of names/header, long sequences to avoid using for merging/gap-filling scaffolds (optional)
 -b  Base name for your output files (optional)
@@ -117,7 +118,7 @@ Usage: ./RAILS [v1.2]
 -f  Assembled Sequences to further scaffold (Multi-FASTA format NO LINE BREAKS, required)
 -q  Long Sequences queried (Multi-FASTA format NO LINE BREAKS, required)
 -s  BAM file (use v1.1 for reading SAM files)
--d  Anchoring bases on contig edges (ie. minimum required alignment size on contigs, default -d 1000, optional)
+-d  Anchoring bases on contig edges (ie. minimum required alignment size on contigs, default -d 250, optional)
 -i  Minimum sequence identity, default -i 0.9, optional
 -t  LIST of names/header, long sequences to avoid using for merging/gap-filling scaffolds (optional)
 -b  Base name for your output files (optional)
@@ -127,7 +128,7 @@ Usage: ./RAILS [v1.2]
 ### How it works
 -------------
 
-The pipeline is detailed in the provided script runRAILS.sh. PLEASE ensure the draft assembly is FASTA-formatted with one sequence per line (NO LINE BREAKS)
+The pipeline is detailed in the provided script runRAILS.sh and runRAILS_lowcontiguityseqs.sh. PLEASE ensure the draft assembly is FASTA-formatted with one sequence per line (NO LINE BREAKS)
 
 Cobbler's process:
 
@@ -136,6 +137,10 @@ In the runRAILS.sh, these scaftigs are renamed, tracking their scaffold of origi
 A bwa index is created and the long sequence file, also re-numbered, is aligned to the scaftigs.
 Cobbler is supplied with the alignment file (-s sam file) and the long reads files (-q option), specifying the minimum length of anchoring bases (-d) aligning at the edge of scaftigs and the minimum sequence identity of the alignment (-i). When 1 or more long sequences align unambiguously to the 3'end of a scaftig and the 5'end of its neighbour, the gap is patched with the sequence of that long sequence. If no long sequences are suitable, or the -d and -i conditions are not met, the original Ns are placed back between those scaftigs.
 
+runRAILS.sh uses scaftigs for patching gaps, whereas runRAILS_lowcontiguityseqs.sh uses scaffold sequences (not broken at Ns) instead.
+If you intend to use an assembly with low contiguity for patching gaps, one that has few and short gaps (stretches of Ns), you may use runRAILS_lowcontiguityseqs.sh.  This is because, depending on the -d parameter set, insufficient anchoring bases may align to patch a gap.
+
+
 RAILS process:
 
 In RAILS, the process is similar as for Cobbler, except that the draft assembly is not broken up at Ns, since the goal is to merge distinct sequences into larger ones.  Long sequences are aligned to the draft assembly sequences, orienting and ordering sequences and simulateneously filling the gaps between them, using DNA bases from the long sequences.

RAILS_v1.2/runRAILS_lowcontiguityseqs.sh

@@ -0,0 +1,31 @@
+#!/bin/bash
+#RLW 2016
+if [ $# -ne 4 ]; then
+        echo "Usage: $(basename $0) <FASTA assembly .fa> <FASTA long sequences .fa> <anchoring sequence length eg. 250> <min sequence identity 0.95>"
+        exit 1
+fi
+###Change line below to point to path of bwa executables
+export PATH=/gsc/btl/linuxbrew/bin:$PATH
+echo Resolving ambiguous bases -Ns- in $1 assembly using long sequences $2
+echo reformatting file $1
+cat $1 | perl -ne 'if(/^\>/){$scafnum++;}else{my $len=length($_);my @scaftigs=split(/N+/i,$_);my $scaftignum=0;foreach my $scaftig(@scaftigs){ my $len=length($scaftig);$scaftignum++;  print ">wga$scafnum";print "."; print "$scaftignum,$len\n$scaftig\n";}}' > $1-formatted.fa
+echo reformatting file $2
+###THIS IS FOR PATCHING GAPS USING SEQUENCES FROM AN ASSEMBLY WITH LOW CONTIGUITY (AND FEW, VERY SHORT STRETCHES OF Ns)
+cat $2 | perl -ne 'if(/^\>/){$ct++;}else{my $len=length($_);print ">seq$ct,$len\n$_";}' > $2-formatted.fa
+echo Building sequence database index out of your $1-formatted.fa assembly contigs..
+bwa index $1-formatted.fa
+echo Aligning long sequences $2-formatted.fa to your contigs..
+bwa mem -a -t24 /projects/btl/rwarren/beluga/RAILS/$1-formatted.fa $2-formatted.fa  | /gsc/btl/linuxbrew/bin/samtools view -Sb - > $2_vs_$1_gapfilling.bam
+echo Scaffolding $1-formatted.fa using $2-formatted.fa and filling gaps with sequences in $2-formatted.fa
+./cobbler.pl -f $1 -s $2_vs_$1_gapfilling.bam -d $3 -i $4 -b $2_vs_$1_$3_gapsFill -q $2-formatted.fa
+echo Process terminated.
+echo RAILS scaffolding $1.gapsFill.fa sequences using long seqs $2 -- anchoring sequence threshold $3 bp 
+echo reformatting file $1.gapsFill.fa
+cat $2_vs_$1_$3_gapsFill.fa | perl -ne 'if(/^\>/){$ct++;}else{my $len=length($_);print ">wga$ct,$len\n$_";}' > $2_vs_$1_$3_gapsFill-formatted.fa
+echo Building sequence database index out of your $2_vs_$1_$3_gapsFill-formatted.fa assembly contigs..
+bwa index $2_vs_$1_$3_gapsFill-formatted.fa
+echo Aligning long sequences $2-formatted.fa to your contigs..
+bwa mem -a -t24 $2_vs_$1_$3_gapsFill-formatted.fa $2-formatted.fa | /gsc/btl/linuxbrew/bin/samtools view -Sb - > $2_vs_$1_scaffolding.bam
+echo Scaffolding $2_vs_$1_$3_gapsFill-formatted.fa using $2-formatted.fa and filling new gaps with sequences in $2-formatted.fa
+./RAILS -f $2_vs_$1_$3_gapsFill-formatted.fa -s $2_vs_$1_scaffolding.bam -d $3 -i $4 -b $2_vs_$1_$3_rails -q $2-formatted.fa &
+echo RAILS process terminated.

rails_v1-2.tar.gz

@@ -6,6 +6,7 @@
 -------------
 
 RAILS: Radial Assembly Improvement by Long Sequence Scaffolding
+
 Cobbler: Gap-filling with long sequences
 
 
@@ -89,7 +90,7 @@ Software. doi: 10.21105/joss.00116
 -------------
 
 <pre>
-./runRAILS.sh
+./runRAILS.sh or runRAILS_lowcontiguityseqs.sh
 Usage: runRAILS.sh <FASTA assembly .fa> <FASTA long sequences .fa> <anchoring sequence length eg. 250> <min sequence identity 0.95>
 
 this pipeline will:
@@ -107,7 +108,7 @@ Usage: ./cobbler.pl [v0.3]
 -f  Assembled Sequences to further scaffold (Multi-FASTA format NO LINE BREAKS, required)
 -q  Long Sequences queried (Multi-FASTA format NO LINE BREAKS, required)
 -s  BAM file (use v0.2 for reading SAM files)
--d  Anchoring bases on contig edges (ie. minimum required alignment size on contigs, default -d 1000, optional)
+-d  Anchoring bases on contig edges (ie. minimum required alignment size on contigs, default -d 250, optional)
 -i  Minimum sequence identity, default -i 0.9, optional
 -t  LIST of names/header, long sequences to avoid using for merging/gap-filling scaffolds (optional)
 -b  Base name for your output files (optional)
@@ -117,7 +118,7 @@ Usage: ./RAILS [v1.2]
 -f  Assembled Sequences to further scaffold (Multi-FASTA format NO LINE BREAKS, required)
 -q  Long Sequences queried (Multi-FASTA format NO LINE BREAKS, required)
 -s  BAM file (use v1.1 for reading SAM files)
--d  Anchoring bases on contig edges (ie. minimum required alignment size on contigs, default -d 1000, optional)
+-d  Anchoring bases on contig edges (ie. minimum required alignment size on contigs, default -d 250, optional)
 -i  Minimum sequence identity, default -i 0.9, optional
 -t  LIST of names/header, long sequences to avoid using for merging/gap-filling scaffolds (optional)
 -b  Base name for your output files (optional)
@@ -127,7 +128,7 @@ Usage: ./RAILS [v1.2]
 ### How it works
 -------------
 
-The pipeline is detailed in the provided script runRAILS.sh. PLEASE ensure the draft assembly is FASTA-formatted with one sequence per line (NO LINE BREAKS)
+The pipeline is detailed in the provided script runRAILS.sh and runRAILS_lowcontiguityseqs.sh. PLEASE ensure the draft assembly is FASTA-formatted with one sequence per line (NO LINE BREAKS)
 
 Cobbler's process:
 
@@ -136,6 +137,10 @@ In the runRAILS.sh, these scaftigs are renamed, tracking their scaffold of origi
 A bwa index is created and the long sequence file, also re-numbered, is aligned to the scaftigs.
 Cobbler is supplied with the alignment file (-s sam file) and the long reads files (-q option), specifying the minimum length of anchoring bases (-d) aligning at the edge of scaftigs and the minimum sequence identity of the alignment (-i). When 1 or more long sequences align unambiguously to the 3'end of a scaftig and the 5'end of its neighbour, the gap is patched with the sequence of that long sequence. If no long sequences are suitable, or the -d and -i conditions are not met, the original Ns are placed back between those scaftigs.
 
+runRAILS.sh uses scaftigs for patching gaps, whereas runRAILS_lowcontiguityseqs.sh uses scaffold sequences (not broken at Ns) instead.
+If you intend to use an assembly with low contiguity for patching gaps, one that has few and short gaps (stretches of Ns), you may use runRAILS_lowcontiguityseqs.sh.  This is because, depending on the -d parameter set, insufficient anchoring bases may align to patch a gap.
+
+
 RAILS process:
 
 In RAILS, the process is similar as for Cobbler, except that the draft assembly is not broken up at Ns, since the goal is to merge distinct sequences into larger ones.  Long sequences are aligned to the draft assembly sequences, orienting and ordering sequences and simulateneously filling the gaps between them, using DNA bases from the long sequences.