Updated WFA2-lib

src/common/wflign/deps/WFA/Makefile

@@ -36,11 +36,11 @@ SUBDIRS=alignment \
         system \
         utils \
         wavefront
-ifeq ($(BUILD_TOOLS),1)
+ifeq ($(BUILD_TOOLS),1)        
     APPS+=tools/generate_dataset \
           tools/align_benchmark
 endif
-ifeq ($(BUILD_EXAMPLES),1)
+ifeq ($(BUILD_EXAMPLES),1)        
     APPS+=examples
 endif
 
@@ -61,13 +61,13 @@ asan: build
 # Build rules
 ###############################################################################
 build: setup
-build: $(SUBDIRS)
-build: lib_wfa
+build: $(SUBDIRS) 
+build: lib_wfa 
 build: $(APPS)
 
 setup:
 	@mkdir -p $(FOLDER_BIN) $(FOLDER_BUILD) $(FOLDER_BUILD_CPP) $(FOLDER_LIB)
-
+	
 lib_wfa: $(SUBDIRS)
 	$(AR) $(AR_FLAGS) $(LIB_WFA) $(FOLDER_BUILD)/*.o 2> /dev/null
 	$(AR) $(AR_FLAGS) $(LIB_WFA_CPP) $(FOLDER_BUILD)/*.o $(FOLDER_BUILD_CPP)/*.o 2> /dev/null
@@ -76,15 +76,16 @@ clean:
 	rm -rf $(FOLDER_BIN) $(FOLDER_BUILD) $(FOLDER_LIB)
 	$(MAKE) --directory=tools/align_benchmark clean
 	$(MAKE) --directory=examples clean
-
+	
 ###############################################################################
 # Subdir rule
 ###############################################################################
 export
 $(SUBDIRS):
 	$(MAKE) --directory=$@ all
-
+	
 $(APPS):
 	$(MAKE) --directory=$@ all
 
 .PHONY: $(SUBDIRS) $(APPS)
+

src/common/wflign/deps/WFA/README.md

@@ -4,7 +4,7 @@
 
 ### 1.1 What is WFA?
 
-The wavefront alignment (WFA) algorithm is an **exact** gap-affine algorithm that takes advantage of homologous regions between the sequences to accelerate the alignment process. Unlike to traditional dynamic programming algorithms that run in quadratic time, the WFA runs in time `O(ns+s^2)`, proportional to the sequence length `n` and the alignment score `s`, using `O(s^2)` memory. Moreover, the WFA algorithm exhibits simple computational patterns that the modern compilers can automatically vectorize for different architectures without adapting the code. To intuitively illustrate why the WFA algorithm is so interesting, look at the following figure. The left panel shows the cells computed by a classical dynamic programming based algorithm (like Smith-Waterman or Needleman Wunsch). In contrast, the right panel shows the cells computed by the WFA algorithm to obtain the same result (i.e., the optimal alignment).
+The wavefront alignment (WFA) algorithm is an **exact** gap-affine algorithm that takes advantage of homologous regions between the sequences to accelerate the alignment process. Unlike to traditional dynamic programming algorithms that run in quadratic time, the WFA runs in time `O(ns+s^2)`, proportional to the sequence length `n` and the alignment score `s`, using `O(s^2)` memory (or `O(s)` using the ultralow/BiWFA mode). Moreover, the WFA algorithm exhibits simple computational patterns that the modern compilers can automatically vectorize for different architectures without adapting the code. To intuitively illustrate why the WFA algorithm is so interesting, look at the following figure. The left panel shows the cells computed by a classical dynamic programming based algorithm (like Smith-Waterman or Needleman Wunsch). In contrast, the right panel shows the cells computed by the WFA algorithm to obtain the same result (i.e., the optimal alignment).
 
 <p align = "center">
 <img src = "img/wfa.vs.swg.png" width="750px">
@@ -37,7 +37,7 @@ Section [WFA2-lib features](#wfa2.features) explores the most relevant options a
     * [Alignment Span](#wfa2.span)
     * [Memory modes](#wfa2.mem)
     * [Heuristic modes](#wfa2.heuristics)
-    <!-- * [Other features](#wfa2.other) -->
+    * [Technical notes](#wfa2.other.notes)
 * [Reporting Bugs and Feature Request](#wfa2.complains)
 * [License](#wfa2.licence)
 * [Citation](#wfa2.cite) 
@@ -47,7 +47,7 @@ Section [WFA2-lib features](#wfa2.features) explores the most relevant options a
 - The WFA algorithm is an **exact algorithm**. If no heuristic is applied (e.g., band or adaptive pruning), the core algorithm guarantees to always find the optimal solution (i.e., best alignment score). Since its first release, some authors have referenced the WFA as approximated or heuristic, which is NOT the case.
 
 
-- Given two sequences of length `n`, traditional dynamic-programming (DP) based methods (like Smith-Waterman or Needleman-Wunsch) compute the optimal alignment in `O(n^2)` time, using `O(n^2)` memory. In contrast, the WFA algorithm requires `O(ns+s^2)` time and `O(s^2)` memory (being `s` the optimal alignment score). Therefore, **the memory consumption of the WFA algorithm is not intrinsically higher than that of other methods**. Most DP-based methods can use heuristics (like banded, X-drop, or Z-drop) to reduce the execution time and the memory usage at the expense of losing accuracy. Likewise, **the WFA algorithm can also use heuristics to reduce the execution time and memory usage**.
+- Given two sequences of length `n`, traditional dynamic-programming (DP) based methods (like Smith-Waterman or Needleman-Wunsch) compute the optimal alignment in `O(n^2)` time, using `O(n^2)` memory. In contrast, the WFA algorithm requires `O(ns+s^2)` time and `O(s^2)` memory (being `s` the optimal alignment score). Therefore, **the memory consumption of the WFA algorithm is not intrinsically higher than that of other methods**. Most DP-based methods can use heuristics (like banded, X-drop, or Z-drop) to reduce the execution time and the memory usage at the expense of losing accuracy. Likewise, **the WFA algorithm can also use heuristics to reduce the execution time and memory usage**. Moreover, the memory mode `ultralow` uses the BiWFA algorithm to execute in `O(ns+s^2)` time and linear `O(s)` memory.
 
 
 - **A note for the fierce competitors.** I can understand that science and publishing have become a fierce competition these days. Many researchers want their methods to be successful and popular, seeking funding, tenure, or even fame. If you are going to benchmark the WFA using the least favourable configuration, careless programming, and a disadvantageous setup, please, go ahead. But remember, researchers like you have put a lot of effort into developing the WFA. We all joined this "competition" because we sought to find better methods that could truly help other researchers. So, try to be nice, tone down the marketing, and produce fair evaluations and honest publications.
@@ -147,7 +147,7 @@ cout << "Alignment score " << aligner.getAlignmentScore() << endl;
 * **Exact alignment** method that computes the optimal **alignment score** and/or **alignment CIGAR**.
 * Supports **multiple distance metrics** (i.e., indel, edit, gap-lineal, gap-affine, and dual-gap gap-affine).
 * Allows performing **end-to-end** (a.k.a. global) and **ends-free** (e.g., semi-global, extension, overlap) alignment.
-* Implements **low-memory modes** to reduce and control memory consumption.
+* Implements **low-memory modes** to reduce and control memory consumption (down to `O(s)` using the `ultralow` mode).
 * Supports various **heuristic strategies** to use on top of the core WFA algorithm.
 * WFA2-lib **operates with plain ASCII strings**. Although we mainly focus on aligning DNA/RNA sequences, the WFA algorithm and the WFA2-lib implementation work with any pair of strings. Moreover, these sequences do not have to be pre-processed (e.g., packed or profiled), nor any table must be precomputed (like the query profile, used within some Smith-Waterman implementations).
 * Due to its simplicity, the WFA algorithm can be automatically vectorized for any SIMD-compliant CPU supported by the compiler. For this reason, **the WFA2-lib implementation is independent of any specific ISA or processor model**. Unlike other hardware-dependent libraries, we aim to offer a multiplatform pairwise-alignment library that can be executed on different processors and models (e.g., SSE, AVX2, AVX512, POWER-ISA, ARM, NEON, SVE, SVE2, RISCV-RVV, ...).
@@ -385,7 +385,7 @@ The WFA2 library allows computing alignments with different spans or shapes. Alt
 
 ### <a name="wfa2.mem"></a> 3.4 Memory modes
 
-The WFA2 library implements various memory modes: `wavefront_memory_high`, `wavefront_memory_med`, `wavefront_memory_low`. These modes allow regulating the overall memory consumption at the expense of execution time. The standard WFA algorithm, which stores explicitly all wavefronts in memory, correspond to the mode `wavefront_memory_high`. The other methods progressively reduce memory usage at the expense of slightly larger alignment times. These memory modes can be used transparently with other alignment options and generate identical results. Note that this option does not affect the score-only alignment mode (it already uses a minimal memory footprint of `O(s)`).
+The WFA2 library implements various memory modes: `wavefront_memory_high`, `wavefront_memory_med`, `wavefront_memory_low`, and , `wavefront_memory_ultralow`. These modes allow regulating the overall memory consumption at the expense of execution time. The standard WFA algorithm, which stores explicitly all wavefronts in memory, correspond to the mode `wavefront_memory_high`. The other methods progressively reduce memory usage at the expense of slightly larger alignment times. These memory modes can be used transparently with other alignment options and generate identical results. Note that this option does not affect the score-only alignment mode (it already uses a minimal memory footprint of `O(s)`).
 
 ```C
   wavefront_aligner_attr_t attributes = wavefront_aligner_attr_default;
@@ -394,7 +394,7 @@ The WFA2 library implements various memory modes: `wavefront_memory_high`, `wave
 
 ### <a name="wfa2.heuristics"></a> 3.5 Heuristic modes
 
-The WFA algorithm can be used combined with many heuristics to reduce the alignment time and memory used. As it happens to other alignment methods, heuristics can result in suboptimal solutions and loss of accuracy. Moreover, some heuristics may drop the alignment if the sequences exceed certain divergence thresholds (i.e., x-drop/z-drop). Due to the popularity and efficiency of these methods, the WFA2 library implements many of these heuristics. Note, **it is not about how little DP-matrix you compute, but about how good the resulting alignments are.**
+The WFA algorithm can be used combined with many heuristics to reduce the alignment time and memory used. As it happens to other alignment methods, heuristics can result in suboptimal solutions and loss of accuracy. Moreover, some heuristics may drop the alignment if the sequences exceed certain divergence thresholds (i.e., x-drop/z-drop). Due to the popularity and efficiency of these methods, the WFA2 library implements many of these heuristics. Note, **it is not about how little DP-matrix you compute, but about how good/accurate the resulting alignments are.**
 
 - **None (for comparison)**. If no heuristic is used, the WFA behaves exploring cells of the DP-matrix in increasing score order (increasing scores correspond to colours from blue to red). 
 
@@ -500,11 +500,20 @@ The WFA algorithm can be used combined with many heuristics to reduce the alignm
   attributes.heuristic.steps_between_cutoffs = 100;
 ```
 
-<!-- ### <a name="wfa2.other"></a> 3.6 Other features -->
+### <a name="wfa2.other.notes"></a> 3.6 Some technical notes
+
+- Thanks to Eizenga's formulation, WFA2-lib can operate with any match score. Although, in practice, M=0 is still the most efficient choice.
+
+- All WFA2-lib algorithms/variants are stable. That is, for alignments having the same score, the library always resolves ties (between M, X, I,and D) using the same criteria: M (highest prio) > X > D > I (lowest prio). Nevertheless, the memory mode `ultralow` (BiWFA) is optimal (always reports the best alignment) but not stable.
 
 ## <a name="wfa2.complains"></a> 4. REPORTING BUGS AND FEATURE REQUEST
 
-Feedback and bug reporting is highly appreciated. Please report any issue or suggestion on github or email to the main developer ([email protected]).
+Feedback and bug reporting is highly appreciated. Please report any issue or suggestion on github or email to the main developer ([email protected]). Don't hesitate to contact us
+if:
+  - You experience any bug or crash.
+  - You want to request a feature or have any suggestion.
+  - Your application using the library is running slower than it should or you expected.
+  - Need help integrating the library into your tool.
 
 ## <a name="wfa2.licence"></a> 5. LICENSE
 
@@ -528,3 +537,5 @@ Miquel Moretó has contributed with fruitful technical discussions and tireless
 
 **Santiago Marco-Sola, Juan Carlos Moure, Miquel Moreto, Antonio Espinosa**. ["Fast gap-affine pairwise alignment using the wavefront algorithm."](https://doi.org/10.1093/bioinformatics/btaa777) Bioinformatics, 2020.
 
+**Santiago Marco-Sola, Jordan M Eizenga, Andrea Guarracino, Benedict Paten, Erik Garrison, Miquel Moreto**. Optimal gap-affine alignment in O(s) space.  _bioRxiv_  (2022). DOI [2022.04.14.488380](https://doi.org/10.1101/2022.04.14.488380)
+

src/common/wflign/deps/WFA/bindings/cpp/WFAligner.cpp

@@ -55,6 +55,7 @@ WFAligner::WFAligner(
     default: this->attributes.memory_mode = wavefront_memory_high; break;
   }
   this->attributes.alignment_scope = (alignmentScope==Score) ? compute_score : compute_alignment;
+  //this->attributes.system.verbose = 2;
   this->wfAligner = nullptr;
 }
 WFAligner::~WFAligner() {
@@ -270,6 +271,19 @@ WFAlignerGapLinear::WFAlignerGapLinear(
   attributes.linear_penalties.indel = indel;
   wfAligner = wavefront_aligner_new(&attributes);
 }
+WFAlignerGapLinear::WFAlignerGapLinear(
+    const int match,
+    const int mismatch,
+    const int indel,
+    const AlignmentScope alignmentScope,
+    const MemoryModel memoryModel) :
+        WFAligner(alignmentScope,memoryModel) {
+  attributes.distance_metric = gap_linear;
+  attributes.linear_penalties.match = match;
+  attributes.linear_penalties.mismatch = mismatch;
+  attributes.linear_penalties.indel = indel;
+  wfAligner = wavefront_aligner_new(&attributes);
+}
 /*
  * Gap-Affine Aligner (a.k.a Smith-Waterman-Gotoh)
  */
@@ -287,6 +301,21 @@ WFAlignerGapAffine::WFAlignerGapAffine(
   attributes.affine_penalties.gap_extension = gapExtension;
   wfAligner = wavefront_aligner_new(&attributes);
 }
+WFAlignerGapAffine::WFAlignerGapAffine(
+    const int match,
+    const int mismatch,
+    const int gapOpening,
+    const int gapExtension,
+    const AlignmentScope alignmentScope,
+    const MemoryModel memoryModel) :
+        WFAligner(alignmentScope,memoryModel) {
+  attributes.distance_metric = gap_affine;
+  attributes.affine_penalties.match = match;
+  attributes.affine_penalties.mismatch = mismatch;
+  attributes.affine_penalties.gap_opening = gapOpening;
+  attributes.affine_penalties.gap_extension = gapExtension;
+  wfAligner = wavefront_aligner_new(&attributes);
+}
 /*
  * Gap-Affine Dual-Cost Aligner (a.k.a. concave 2-pieces)
  */
@@ -308,5 +337,25 @@ WFAlignerGapAffine2Pieces::WFAlignerGapAffine2Pieces(
   attributes.affine2p_penalties.gap_extension2 = gapExtension2;
   wfAligner = wavefront_aligner_new(&attributes);
 }
+WFAlignerGapAffine2Pieces::WFAlignerGapAffine2Pieces(
+    const int match,
+    const int mismatch,
+    const int gapOpening1,
+    const int gapExtension1,
+    const int gapOpening2,
+    const int gapExtension2,
+    const AlignmentScope alignmentScope,
+    const MemoryModel memoryModel) :
+        WFAligner(alignmentScope,memoryModel) {
+  attributes.distance_metric = gap_affine_2p;
+  attributes.affine2p_penalties.match = match;
+  attributes.affine2p_penalties.mismatch = mismatch;
+  attributes.affine2p_penalties.gap_opening1 = gapOpening1;
+  attributes.affine2p_penalties.gap_extension1 = gapExtension1;
+  attributes.affine2p_penalties.gap_opening2 = gapOpening2;
+  attributes.affine2p_penalties.gap_extension2 = gapExtension2;
+  wfAligner = wavefront_aligner_new(&attributes);
+}
 
 } /* namespace wfa */
+

src/common/wflign/deps/WFA/bindings/cpp/WFAligner.hpp

@@ -61,7 +61,7 @@ class WFAligner {
   };
   enum AlignmentStatus {
     StatusSuccessful = WF_STATUS_SUCCESSFUL,
-    StatusDropped = WF_STATUS_HEURISTICALY_DROPPED,
+    StatusUnfeasible = WF_STATUS_UNFEASIBLE,
     StatusMaxScoreReached = WF_STATUS_MAX_SCORE_REACHED,
     StatusOOM = WF_STATUS_OOM,
   };
@@ -183,6 +183,12 @@ class WFAlignerGapLinear : public WFAligner {
       const int indel,
       const AlignmentScope alignmentScope,
       const MemoryModel memoryModel = MemoryHigh);
+  WFAlignerGapLinear(
+      const int match,
+      const int mismatch,
+      const int indel,
+      const AlignmentScope alignmentScope,
+      const MemoryModel memoryModel = MemoryHigh);
 };
 /*
  * Gap-Affine Aligner (a.k.a Smith-Waterman-Gotoh)
@@ -195,6 +201,13 @@ class WFAlignerGapAffine : public WFAligner {
       const int gapExtension,
       const AlignmentScope alignmentScope,
       const MemoryModel memoryModel = MemoryHigh);
+  WFAlignerGapAffine(
+      const int match,
+      const int mismatch,
+      const int gapOpening,
+      const int gapExtension,
+      const AlignmentScope alignmentScope,
+      const MemoryModel memoryModel = MemoryHigh);
 };
 /*
  * Gap-Affine Dual-Cost Aligner (a.k.a. concave 2-pieces)
@@ -209,6 +222,15 @@ class WFAlignerGapAffine2Pieces : public WFAligner {
       const int gapExtension2,
       const AlignmentScope alignmentScope,
       const MemoryModel memoryModel = MemoryHigh);
+  WFAlignerGapAffine2Pieces(
+      const int match,
+      const int mismatch,
+      const int gapOpening1,
+      const int gapExtension1,
+      const int gapOpening2,
+      const int gapExtension2,
+      const AlignmentScope alignmentScope,
+      const MemoryModel memoryModel = MemoryHigh);
 };
 
 } /* namespace wfa */

src/common/wflign/deps/WFA/tools/align_benchmark/align_benchmark.c

@@ -139,7 +139,7 @@ typedef struct {
 } benchmark_args;
 benchmark_args parameters = {
   // Algorithm
-  .algorithm = alignment_test,
+  .algorithm = alignment_edit_wavefront,
   // I/O
   .input_filename = NULL,
   .output_filename = NULL,
@@ -1060,8 +1060,8 @@ void parse_arguments(int argc,char** argv) {
         parameters.verbose = 1;
       } else {
         parameters.verbose = atoi(optarg);
-        if (parameters.verbose < 0 || parameters.verbose > 3) {
-          fprintf(stderr,"Option '--verbose' must be in {0,1,2,3}\n");
+        if (parameters.verbose < 0 || parameters.verbose > 4) {
+          fprintf(stderr,"Option '--verbose' must be in {0,1,2,3,4}\n");
           exit(1);
         }
       }