Allow skipping out of loop with LOCO HEAD master

2 files changed, 21 insertions, 27 deletions

diff --git a/src/lmm.cpp b/src/lmm.cpp
index 969e6fc..2b730f3 100644
--- a/src/lmm.cpp
+++ b/src/lmm.cpp
@@ -2012,23 +2012,17 @@ void LMM::mdb_analyze(std::function< SnpNameValues2(size_t) >& fetch_snp,
     // continue;
 
     auto tup = fetch_snp(t); // use the callback
-    auto success = get<0>(tup);
-    if (!success)
+    auto state = get<0>(tup);
+    if (state == SKIP)
       continue;
-    // typedef tuple< bool,MarkerInfo,vector<double> > SnpNameValues2;
-    // auto marker = get<1>(tup);
-    // auto chr = get<2>(tup);
-    // auto mpos = get<3>(tup);
+    if (state == LAST)
+      break; // marker loop because of LOCO
+
     auto markerinfo = get<1>(tup);
     auto gs = get<2>(tup);
 
     markers.push_back(markerinfo);
 
-    // check whether SNP is included in gwasnps (used by LOCO)
-    /*
-    if (process_gwasnps && gwasnps.count(snp) == 0)
-      continue;
-    */
     // drop missing idv and plug mean values for missing geno
     double x_total = 0.0; // sum genotype values to compute x_mean
     uint vpos = 0;        // position in target vector
@@ -2067,7 +2061,6 @@ void LMM::mdb_analyze(std::function< SnpNameValues2(size_t) >& fetch_snp,
     gsl_vector_safe_memcpy(&Xlarge_col.vector, x);
     c++; // count markers going in
 
-
     if (c % msize == 0) {
       batch_compute(msize,markers);
       markers.clear();
@@ -2212,10 +2205,6 @@ void LMM::mdb_calc_gwa(const gsl_matrix *U, const gsl_vector *eval,
   mdb_stat(rtxn, geno_mdb, &stat);
   auto num_markers = stat.ms_entries;
 
-  // fetch_snp is a callback function for every SNP row
-  // returns typedef std::tuple<bool,string,std::vector<double> > SnpNameValues;
-
-  // size_t prev_line = 0;
   auto mdb_fetch = MDB_FIRST;
 
   auto cursor = lmdb::cursor::open(rtxn, geno_mdb);
@@ -2227,20 +2216,14 @@ void LMM::mdb_calc_gwa(const gsl_matrix *U, const gsl_vector *eval,
 
     string_view key,value;
 
-    /*
-    while (prev_line <= num) {
-      cursor.get(key, value, MDB_NEXT);
-      prev_line++;
-    }
-    */
-    auto success = cursor.get(key, value, mdb_fetch);
+    auto mdb_success = cursor.get(key, value, mdb_fetch);
     mdb_fetch = MDB_NEXT;
 
     // uint8_t chr;
     vector<double> gs;
     MarkerInfo markerinfo;
 
-    if (success) {
+    if (mdb_success) {
       size_t size = 0;
       // ---- Depending on the format we get different buffers - currently float and byte are supported:
       if (format == "Gb") {
@@ -2282,7 +2265,10 @@ void LMM::mdb_calc_gwa(const gsl_matrix *U, const gsl_vector *eval,
       // printf("%#02x %#02x\n", chr, loco_chr);
 
       if (is_loco && loco_chr != chr) {
-        return make_tuple(false, MarkerInfo { .name="", .chr=chr, .pos=pos } , gs);
+        if (chr > loco_chr)
+            return make_tuple(LAST, MarkerInfo { .name="", .chr=chr, .pos=pos } , gs);
+          else
+            return make_tuple(SKIP, MarkerInfo { .name="", .chr=chr, .pos=pos } , gs);
       }
 
       string_view value2;
@@ -2299,7 +2285,7 @@ void LMM::mdb_calc_gwa(const gsl_matrix *U, const gsl_vector *eval,
 
       // cout << "!!!!" << size << marker << ": af" << maf << " " << gs[0] << "," << gs[1] << "," << gs[2] << "," << gs[3] << endl;
     }
-    return make_tuple(success, markerinfo, gs);
+    return make_tuple(COMPUTE, markerinfo, gs);
   };
   LMM::mdb_analyze(fetch_snp,U,eval,UtW,Uty,W,y,num_markers);
 
diff --git a/src/lmm.h b/src/lmm.h
index da5ad21..295602a 100644
--- a/src/lmm.h
+++ b/src/lmm.h
@@ -55,7 +55,15 @@ struct MarkerInfo {
 
 typedef vector<MarkerInfo> Markers;
 typedef tuple< string,vector<double> > SnpNameValues;
-typedef tuple< bool,MarkerInfo,vector<double> > SnpNameValues2; // success, markerinfo (maf and n_miss are computed)
+
+enum MarkerState {
+  FAIL,
+  COMPUTE,
+  SKIP,
+  LAST
+};
+
+typedef tuple< MarkerState,MarkerInfo,vector<double> > SnpNameValues2; // success, markerinfo (maf and n_miss are computed)
 // Results for LMM.
 class SUMSTAT2 {
 public: