Added tabular reporting option.

Author: bo1929

src/common.hpp

@@ -45,7 +45,7 @@
   #include <immintrin.h>
 #endif
 
-#define STRSTREAM_PRECISION 6
+#define STRSTREAM_PRECISION 4
 
 #define VERSION "v0.6.0"
 #define PRINT_VERSION std::cerr << "krepp version: " << VERSION << std::endl;

src/krepp.cpp

@@ -403,6 +403,8 @@ void QueryIndex::header_preport(strstream& dreport_stream)
   qtree->stream_nwk_jplace(dreport_stream, qtree->get_root());
   if (summarize) {
     dreport_stream << "\nDISTAL_NODE\tEDGE_NUM\tWEIGHTED_COUNT\tSEQUENCE_ABUNDANCE\n";
+  } else if (tabular) {
+    dreport_stream << "\nSEQ_ID\tDISTAL_NODE\tEDGE_NUM\tLWR\tDIST\n";
   } else {
     (void)0; // TODO: Perhaps introduce an alternative placement report format.
   }
@@ -439,7 +441,7 @@ void QueryIndex::place_sequences()
   strstream preport_stream;
   preport_stream.precision(STRSTREAM_PRECISION);
   preport_stream << std::fixed;
-  if (summarize) {
+  if (summarize || tabular) {
     header_preport(preport_stream);
     (*output_stream) << preport_stream.rdbuf();
   } else {
@@ -461,7 +463,7 @@ void QueryIndex::place_sequences()
 #pragma omp task untied
         {
           strstream batch_stream;
-          ib.place_sequences(batch_stream);
+          ib.place_sequences(batch_stream, tabular);
           if (!summarize && !cont_reading) {
 #pragma omp task
             {
@@ -498,6 +500,8 @@ void QueryIndex::place_sequences()
       preport_stream << nd->get_name(true) << "\t" << nd->get_en() << "\t" << wcount << "\t" << wcount / twcount << "\n";
     }
     (*output_stream) << preport_stream.rdbuf();
+  } else if (tabular) {
+    (void)0;
   } else {
     preport_stream.str("");
     end_jplace(preport_stream);
@@ -615,6 +619,9 @@ void QueryIndex::init_sc_place(CLI::App& sc)
   sc.add_flag("--filter,!--no-filter",
               filter,
               "Filter a placement when there is not enough k-mer matches below threshold tau. [true]");
+  sc.add_flag("--tabular,!--no-tabular",
+              tabular,
+              "Output the per query sequence placements (taxonomic/phylogenetic) in a tab-separated format. [false]");
   sc.callback([&]() {
     no_filter = !filter;
     if (!output_path.empty()) {
@@ -747,12 +754,11 @@ int main(int argc, char** argv)
     if (sc_place.count("--lineage-file") + sc_place.count("-l")) {
       krepp_place.read_lineages();
       std::cerr << "Placing given sequences on the taxonomic lineage..." << std::endl;
-      krepp_place.place_sequences();
     } else {
       krepp_place.ensure_backbone();
       std::cerr << "Placing given sequences on the backbone tree..." << std::endl;
-      krepp_place.place_sequences();
     }
+    krepp_place.place_sequences();
     std::chrono::duration<float> es_s = std::chrono::system_clock::now() - tstart - es_b;
     std::cerr << "Done placing queries, elapsed: " << es_s.count() << " sec" << std::endl;
     std::cerr << "Total number of sequences queried: " << krepp_place.get_total_qseq() << std::endl;

src/krepp.hpp

@@ -212,6 +212,7 @@ class QueryIndex : public TargetIndex
   double dist_max = 0.25;
   bool no_filter = true;
   bool filter = false;
+  bool tabular = false;
   bool multi = true;
   bool summarize = false;
   uint64_t total_qseq = 0;

src/query.cpp

@@ -192,7 +192,7 @@ void IBatch::report_distances(strstream& batch_stream)
   }
 }
 
-void IBatch::place_sequences(strstream& batch_stream)
+void IBatch::place_sequences(strstream& batch_stream, bool tabular)
 {
   for (bix = 0; bix < batch_size; ++bix) {
     const char* seq = seq_batch[bix].data();
@@ -205,11 +205,11 @@ void IBatch::place_sequences(strstream& batch_stream)
     summarize_matches(imers_or, imers_rc);
     batch_stream.precision(STRSTREAM_PRECISION);
     batch_stream << std::fixed;
-    report_placement(batch_stream);
+    report_placement(batch_stream, tabular);
   }
 }
 
-void IBatch::report_placement(strstream& batch_stream)
+void IBatch::report_placement(strstream& batch_stream, bool tabular)
 {
   if (node_to_minfo.size() == 0 || !(no_filter || (mi_closest->get_leq_tau(tau) > 1.0))) {
     return;
@@ -218,12 +218,16 @@ void IBatch::report_placement(strstream& batch_stream)
   minfo_sptr_t mi_pp = mi_closest;
   mi_pp->chisq = 0;
 
-  batch_stream << "\t\t\t{\"n\" : [\"" << identifer_batch[bix] << "\"], \"p\" : [";
+  if (!tabular) batch_stream << "\t\t\t{\"n\" : [\"" << identifer_batch[bix] << "\"], \"p\" : [";
   if (node_to_minfo.size() == 1) {
     if (summarize) {
       node_to_wcount[nd_pp] += 1.0;
     } else {
-      batch_stream << PLACEMENT_FIELD(nd_pp, mi_pp) << "]},\n";
+      if (tabular) {
+        batch_stream << identifer_batch[bix] << "\t" << TABULAR_FIELD(nd_pp, mi_pp) << "\n";
+      } else {
+        batch_stream << PLACEMENT_FIELD(nd_pp, mi_pp) << "]},\n";
+      }
     }
     return;
   }
@@ -281,11 +285,15 @@ void IBatch::report_placement(strstream& batch_stream)
       if (summarize) {
         node_to_wcount[nd_pp] += 1.0 / nd_v.size();
       } else {
-        if (i > 0) batch_stream << ",";
-        batch_stream << "\n\t\t\t\t" << PLACEMENT_FIELD(nd_pp, mi_pp);
+        if (i > 0 && !tabular) batch_stream << ",";
+        if (tabular) {
+          batch_stream << identifer_batch[bix] << "\t" << TABULAR_FIELD(nd_pp, mi_pp) << "\n";
+        } else {
+          batch_stream << "\n\t\t\t\t" << PLACEMENT_FIELD(nd_pp, mi_pp);
+        }
       }
     }
-    if (!summarize) {
+    if (!summarize && !tabular) {
       batch_stream << "]\n\t\t\t},\n";
     }
   } else {
@@ -302,7 +310,11 @@ void IBatch::report_placement(strstream& batch_stream)
     if (summarize) {
       node_to_wcount[nd_pp] += 1.0;
     } else {
-      batch_stream << PLACEMENT_FIELD(nd_pp, mi_pp) << "]},\n";
+      if (tabular) {
+        batch_stream << identifer_batch[bix] << "\t" << TABULAR_FIELD(nd_pp, mi_pp) << "\n";
+      } else {
+        batch_stream << PLACEMENT_FIELD(nd_pp, mi_pp) << "]},\n";
+      }
     }
   }
 }

src/query.hpp

@@ -58,9 +58,9 @@ class IBatch
   void search_mers(const char* seq, uint64_t len, imers_sptr_t imers_or, imers_sptr_t imers_rc);
   void summarize_matches(imers_sptr_t imers_or, imers_sptr_t imers_rc);
   void estimate_distances(strstream& batch_stream);
-  void place_sequences(strstream& batch_stream);
   void report_distances(strstream& batch_stream);
-  void report_placement(strstream& batch_stream);
+  void place_sequences(strstream& batch_stream, bool tabular);
+  void report_placement(strstream& batch_stream, bool tabular);
   const parallel_flat_phmap<node_sptr_t, double>& get_summary() { return node_to_wcount; }
 
 private:
@@ -191,6 +191,8 @@ class Minfo
   "[" << nd->get_en() << ", " << mi->jukes_cantor_dist() - nd->get_midpoint_pendant() << ", " << nd->get_midpoint_pendant() \
       << ", " << -mi->v_llh << ", " << mi->lwr << ", " << mi->d_llh << "]"
 
+#define TABULAR_FIELD(nd, mi) nd->get_name(true) << "\t" << nd->get_en() << "\t" << mi->lwr << "\t" << mi->d_llh
+
 #define DISTANCE_FIELD(nd, mi) nd->get_name() << "\t" << mi->d_llh << "\n"
 
 private: