Count-based differential assembly

metagraph/.clang-format → .clang-format

@@ -41,7 +41,7 @@ ColumnLimit: 90
 CompactNamespaces: false
 ConstructorInitializerAllOnOneLineOrOnePerLine: true
 ConstructorInitializerIndentWidth: 4
-ContinuationIndentWidth: 8
+ContinuationIndentWidth: 4
 Cpp11BracedListStyle: false
 DerivePointerAlignment: false
 DisableFormat: false
@@ -51,11 +51,13 @@ ForEachMacros:
   - foreach
   - Q_FOREACH
   - BOOST_FOREACH
+IncludeBlocks: Regroup
 IndentCaseLabels: true
 IndentPPDirectives: None
 IndentWidth: 4
 IndentWrappedFunctionNames: false
 KeepEmptyLinesAtTheStartOfBlocks: false
+LambdaBodyIndentation: Signature
 MaxEmptyLinesToKeep: 2
 NamespaceIndentation: None
 PenaltyBreakAssignment: 2

.gitmodules

@@ -13,7 +13,7 @@
 	url = https://github.com/karasikov/caches.git
 [submodule "metagraph/external-libraries/eigen"]
 	path = metagraph/external-libraries/eigen
-	url = https://github.com/eigenteam/eigen-git-mirror.git
+	url = https://gitlab.com/libeigen/eigen.git
 [submodule "metagraph/external-libraries/KMC"]
 	path = metagraph/external-libraries/KMC
 	url = https://github.com/karasikov/KMC.git

metagraph/src/annotation/int_matrix/base/int_matrix.cpp

@@ -1,5 +1,11 @@
 #include "int_matrix.hpp"
 
+#include <progress_bar.hpp>
+
+#include "common/logger.hpp"
+#include "common/threads/threading.hpp"
+#include "common/vectors/vector_algorithm.hpp"
+
 
 namespace mtg {
 namespace annot {
@@ -8,6 +14,75 @@ namespace matrix {
 using Row = BinaryMatrix::Row;
 using Column = BinaryMatrix::Column;
 
+void IntMatrix::call_row_values(const std::function<void(uint64_t, const RowValues&, size_t)> &callback,
+                                bool ordered) const {
+    size_t n = get_binary_matrix().num_rows();
+    size_t batch_size = (n + get_num_threads() - 1) / get_num_threads();
+    size_t rows_per_update = 10000;
+    ProgressBar progress_bar(n, "Streaming rows", std::cerr, !common::get_verbose());
+
+    std::ignore = ordered;
+
+    #pragma omp parallel for num_threads(get_num_threads()) schedule(dynamic)
+    for (size_t k = 0; k < n; k += batch_size) {
+        size_t begin = k;
+        size_t end = std::min(begin + batch_size, n);
+        size_t thread_id = k / batch_size;
+        std::vector<Row> row_ids;
+        for (size_t row_batch_i = begin; row_batch_i < end; row_batch_i += rows_per_update) {
+            size_t row_batch_end = std::min(row_batch_i + rows_per_update, end);
+            row_ids.resize(row_batch_end - row_batch_i);
+            std::iota(row_ids.begin(), row_ids.end(), row_batch_i);
+            auto row_vals = get_row_values(row_ids);
+            for (size_t i = 0; i < row_ids.size(); ++i) {
+                callback(i + row_batch_i, row_vals[i], thread_id);
+            }
+            progress_bar += row_ids.size();
+        }
+    }
+}
+
+std::vector<VectorMap<uint64_t, size_t>> IntMatrix::get_histograms(const std::vector<size_t> &min_counts,
+                                                                   sdsl::bit_vector *unmark_discarded) const {
+    common::logger->trace("Calculating count histograms");
+    bool parallel = get_num_threads() > 0;
+    std::vector<VectorMap<uint64_t, size_t>> hists_map(get_binary_matrix().num_columns());
+    std::atomic_thread_fence(std::memory_order_release);
+    call_row_values([&](uint64_t row_i, const auto &row, size_t) {
+        if (min_counts.size()) {
+            bool keep_row = false;
+            for (const auto &[j, c] : row) {
+                if (c >= min_counts[j]) {
+                    keep_row = true;
+                    break;
+                }
+            }
+
+            if (!keep_row) {
+                if (unmark_discarded)
+                    unset_bit(unmark_discarded->data(), row_i, parallel, std::memory_order_relaxed);
+
+                return;
+            }
+        }
+
+        Vector<uint64_t> counts(hists_map.size());
+        for (const auto &[j, raw_c] : row) {
+            if (min_counts.empty() || raw_c >= min_counts[j])
+                counts[j] = raw_c;
+        }
+
+        #pragma omp critical
+        {
+            for (size_t j = 0; j < counts.size(); ++j) {
+                ++hists_map[j][counts[j]];
+            }
+        }
+    }, false);
+    std::atomic_thread_fence(std::memory_order_acquire);
+
+    return hists_map;
+}
 
 IntMatrix::RowValues
 IntMatrix::sum_row_values(const std::vector<std::pair<Row, size_t>> &index_counts,

metagraph/src/annotation/int_matrix/base/int_matrix.hpp

@@ -3,7 +3,10 @@
 
 #include <vector>
 
+#include <sdsl/int_vector.hpp>
+
 #include "annotation/binary_matrix/base/binary_matrix.hpp"
+#include "common/vector_map.hpp"
 
 
 namespace mtg {
@@ -20,6 +23,12 @@ class IntMatrix {
     // |row| is in [0, num_rows), |column| is in [0, num_columns)
     virtual std::vector<RowValues> get_row_values(const std::vector<BinaryMatrix::Row> &rows) const = 0;
 
+    virtual void call_row_values(const std::function<void(uint64_t, const RowValues&, size_t)> &callback,
+                                 bool ordered = true) const;
+
+    virtual std::vector<VectorMap<uint64_t, size_t>> get_histograms(const std::vector<size_t> &min_counts = {},
+                                                                    sdsl::bit_vector *unmark_discarded = nullptr) const;
+
     // sum up values for each column with at least |min_count| non-zero values
     virtual RowValues
     sum_row_values(const std::vector<std::pair<BinaryMatrix::Row, size_t>> &index_counts,