start with eaf get stuck at local optimal

Zilong-Li · Jan 31, 2024 · 0fcf704 · 0fcf704
1 parent 45aba95
commit 0fcf704
Show file tree

Hide file tree

Showing 5 changed files with 2 additions and 19 deletions.
diff --git a/Makefile b/Makefile
@@ -3,7 +3,7 @@ CXX      = g++
 
 # CXXFLAGS = -std=c++17 -Wall -O3 -g -fsanitize=address
 # CXXFLAGS = -std=c++17 -Wall -O3 -march=native -DNDEBUG
-CXXFLAGS = -std=c++17 -Wall -O3  -mavx2 -fPIC -DNDEBUG
+CXXFLAGS = -std=c++17 -Wall -O3  -mavx2
 INC      = -I./src -I./inst/include -I$(HTSDIR)
 LDFLAGS  =  -L$(HTSDIR) -Wl,-rpath,$(HTSDIR)
 LIBS     = $(HTSDIR)/libhts.a -llzma -lbz2 -lm -lz -lpthread

diff --git a/src/common.hpp b/src/common.hpp
@@ -97,7 +97,7 @@ struct Options
     double ftol{1e-6}; // threshold for F
     double qtol{1e-6}; // threshold for Q
     bool noaccel{0}, noscreen{0}, single_chunk{0}, debug{0}, collapse{0};
-    bool nQ{0}, nP{0}, nF{0}, nR{0}, aQ{0}, oVCF{0}, eAF{1}, eHap{1};
+    bool nQ{0}, nP{0}, nF{0}, nR{0}, aQ{0}, oVCF{0}, eHap{0};
     std::string out, in_beagle, in_vcf, in_bin, in_impute, in_joint;
     std::string samples{""}, region{""}, in_plink{""}, in_qfile{""}, in_pfile{""}, in_rfile{""};
     std::string opts_in_effect{"Options in effect:\n   "};

diff --git a/src/fastphase.cpp b/src/fastphase.cpp
@@ -42,16 +42,6 @@ void FastPhaseK2::setFlags(double tol_p, double tol_f, double tol_q, bool debug_
     NR = nR;
 }
 
-void FastPhaseK2::setStartPoint(const std::unique_ptr<BigAss> & genome)
-{
-    for(size_t ic = 0; ic < genome->pos.size(); ic++)
-    {
-        MyArr1D af = estimate_af_by_gl(genome->gls[ic], genome->nsamples, genome->pos[ic].size()).cast<MyFloat>();
-        const int S = pos_chunk[ic + 1] - pos_chunk[ic];
-        for(int s = 0; s < S; s++) F.row(pos_chunk[ic] + s) = af(s);
-    }
-}
-
 void FastPhaseK2::initIteration()
 {
     // initial temp variables
@@ -216,7 +206,6 @@ int run_impute_main(Options & opts)
     FastPhaseK2 faith(genome->nsamples, genome->nsnps, opts.C, opts.seed);
     faith.setFlags(opts.ptol, opts.ftol, opts.qtol, opts.debug, opts.nQ, opts.nP, opts.nF, opts.nR);
     faith.initRecombination(genome->pos, opts.in_rfile);
-    if(opts.eAF) faith.setStartPoint(genome);
     double loglike, diff, prevlike{std::numeric_limits<double>::lowest()};
     for(int it = 0; SIG_COND && it <= opts.nimpute; it++)
     {

diff --git a/src/fastphase.hpp b/src/fastphase.hpp
@@ -55,7 +55,6 @@ class FastPhaseK2
 
     void initRecombination(const Int2D & pos, std::string rfile = "", int B = 1, double Ne = 20000);
     void setFlags(double, double, double, bool, bool, bool, bool, bool);
-    void setStartPoint(const std::unique_ptr<BigAss> & genome);
     void protectPars();
     void initIteration();
     void updateIteration();

diff --git a/src/main.cpp b/src/main.cpp
@@ -176,10 +176,6 @@ int main(int argc, char * argv[])
         .help("seed for reproducibility")
         .default_value(999)
         .scan<'i', int>();
-    cmd_impute.add_argument("--no-eaf")
-        .help("do not start with eaf")
-        .default_value(false)
-        .implicit_value(true);
     cmd_impute.add_argument("--write-hapsum")
         .help("write Hapsum instead of AE into parse.bin")
         .default_value(false)
@@ -295,7 +291,6 @@ int main(int argc, char * argv[])
             opts.seed = cmd_impute.get<int>("--seed");
             opts.chunksize = cmd_impute.get<int>("--chunksize");
             opts.single_chunk = cmd_impute.get<bool>("--single-chunk");
-            opts.eAF = !cmd_impute.get<bool>("--no-eaf");
             opts.eHap = cmd_impute.get<bool>("--write-hapsum");
             opts.collapse = cmd_impute.get<bool>("--collapse");
             opts.tol_r = cmd_impute.get<double>("--minRecombRate");