[Genabel-commits] r889 - branches/ProbABEL-refactoring/ProbABEL/src

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Mon Apr 9 15:27:08 CEST 2012


Author: maartenk
Date: 2012-04-09 15:27:08 +0200 (Mon, 09 Apr 2012)
New Revision: 889

Modified:
   branches/ProbABEL-refactoring/ProbABEL/src/coxph_data.cpp
   branches/ProbABEL-refactoring/ProbABEL/src/data.cpp
   branches/ProbABEL-refactoring/ProbABEL/src/main.cpp
Log:
-added brackets to cmpfun in coxph_data.cpp  
data.cpp:
- removed  ncov and ndis because they were unused variables
- removed extern bool is_interaction_excluded
- made scoop of row col and val smaller in InvSigma::InvSigma 
main.cpp
-removed unused variable null_loglik
-extracted code to separate function:
 open_files_for_output() abstracts code  and code duplication
update_progress_to_cmd_line() abstracts 

Modified: branches/ProbABEL-refactoring/ProbABEL/src/coxph_data.cpp
===================================================================
--- branches/ProbABEL-refactoring/ProbABEL/src/coxph_data.cpp	2012-04-08 14:37:34 UTC (rev 888)
+++ branches/ProbABEL-refactoring/ProbABEL/src/coxph_data.cpp	2012-04-09 13:27:08 UTC (rev 889)
@@ -21,11 +21,17 @@
     double el1 = *(double*) a;
     double el2 = *(double*) b;
     if (el1 > el2)
+    {
         return 1;
+    }
     if (el1 < el2)
+    {
         return -1;
+    }
     if (el1 == el2)
+    {
         return 0;
+    }
 }
 
 coxph_data::coxph_data(const coxph_data &obj)

Modified: branches/ProbABEL-refactoring/ProbABEL/src/data.cpp
===================================================================
--- branches/ProbABEL-refactoring/ProbABEL/src/data.cpp	2012-04-08 14:37:34 UTC (rev 888)
+++ branches/ProbABEL-refactoring/ProbABEL/src/data.cpp	2012-04-09 13:27:08 UTC (rev 889)
@@ -20,14 +20,13 @@
 #include "utilities.h"
 using namespace std;
 
-extern bool is_interaction_excluded;
-
 unsigned int Nmeasured(char * fname, int nphenocols, int npeople)
 {
-    int ncov = nphenocols - 2;
-    int nids_all = npeople;
+//TODO: unused variables remove them for good if there is no reason to keep them
+//int ncov = nphenocols - 2;
+//int nids_all = npeople;
 
-    // first pass -- find unmeasured people
+// first pass -- find unmeasured people
     std::ifstream infile(fname);
     if (!infile)
     {
@@ -170,9 +169,7 @@
     npeople = phe->nids;
     std::ifstream myfile(filename_);
     char * line = new char[MAXIMUM_PEOPLE_AMOUNT];
-    double val;
     std::string id;
-    unsigned row = 0, col = 0;
 
     matrix.reinit(npeople, npeople);
 
@@ -180,6 +177,8 @@
 
     if (myfile.is_open())
     {
+        double val;
+        unsigned row = 0;
         while (myfile.getline(line, MAXIMUM_PEOPLE_AMOUNT))
         {
 
@@ -194,7 +193,7 @@
                         << " must be there. Wrong inverse variance matrix (only measured id must be there)\n";
                 exit(1);
             }
-
+            unsigned col = 0;
             while (line_stream >> val)
             {
                 matrix.put(val, row, col);

Modified: branches/ProbABEL-refactoring/ProbABEL/src/main.cpp
===================================================================
--- branches/ProbABEL-refactoring/ProbABEL/src/main.cpp	2012-04-08 14:37:34 UTC (rev 888)
+++ branches/ProbABEL-refactoring/ProbABEL/src/main.cpp	2012-04-09 13:27:08 UTC (rev 889)
@@ -46,6 +46,70 @@
 #define EPS 1.e-8
 #define CHOLTOL 1.5e-12
 
+void update_progress_to_cmd_line(int csnp, int nsnps)
+{
+    if (csnp % 1000 == 0)
+    {
+        if (csnp == 0)
+        {
+            fprintf(stdout, "Analysis: %6.2f ...",
+                    100. * double(csnp) / double(nsnps));
+        }
+        else
+        {
+            fprintf(stdout, "\b\b\b\b\b\b\b\b\b\b%6.2f ...",
+                    100. * double(csnp) / double(nsnps));
+        }
+        std::cout.flush();
+    }
+}
+
+void open_files_for_output(std::vector<std::ofstream*>& outfile,
+        std::string& outfilename_str)
+{
+    // open a file for output
+    //_____________________
+    for (int i = 0; i < 5; i++)
+    {
+        outfile.push_back(new std::ofstream());
+    }
+    outfile[0]->open((outfilename_str + "_2df.out.txt").c_str());
+    outfile[1]->open((outfilename_str + "_add.out.txt").c_str());
+    outfile[2]->open((outfilename_str + "_domin.out.txt").c_str());
+    outfile[3]->open((outfilename_str + "_recess.out.txt").c_str());
+    outfile[4]->open((outfilename_str + "_over_domin.out.txt").c_str());
+    if (!outfile[0]->is_open())
+    {
+        std::cerr << "Cannot open file for writing: "
+                << outfilename_str + "_2df.out.txt" << "\n";
+        exit(1);
+    }
+    if (!outfile[1]->is_open())
+    {
+        std::cerr << "Cannot open file for writing: "
+                << outfilename_str + "_add.out.txt" << "\n";
+        exit(1);
+    }
+    if (!outfile[2]->is_open())
+    {
+        std::cerr << "Cannot open file for writing: "
+                << outfilename_str + "_domin.out.txt" << "\n";
+        exit(1);
+    }
+    if (!outfile[3]->is_open())
+    {
+        std::cerr << "Cannot open file for writing: "
+                << outfilename_str + "_recess.out.txt" << "\n";
+        exit(1);
+    }
+    if (!outfile[4]->is_open())
+    {
+        std::cerr << "Cannot open file for writing: "
+                << outfilename_str + "_over_domin.out.txt" << "\n";
+        exit(1);
+    }
+}
+
 int main(int argc, char * argv[])
 {
 
@@ -160,11 +224,12 @@
      gendata gtd (input_var.getGenfilename(),nsnps,input_var.getNgpreds(),phd.nids_all,phd.nids,phd.allmeasured,skipd,phd.idnames);
      **/
     // estimate null model
-    double null_loglik = 0.;
+    //TODO: remove this unused variable if there is not a reason to keep it
+    //double null_loglik = 0.;
 #if COXPH
     coxph_data nrgd=coxph_data(phd,gtd,-1,input_var.isIsInteractionExcluded());
 #else
-    regdata nrgd = regdata(phd, gtd, -1,input_var.isIsInteractionExcluded());
+    regdata nrgd = regdata(phd, gtd, -1, input_var.isIsInteractionExcluded());
 #endif
 
     std::cout << " loaded null data ...";
@@ -183,7 +248,7 @@
 
     nrd.estimate(nrgd,0,MAXITER,EPS,CHOLTOL,0, input_var.getInteraction(), input_var.getNgpreds(), 1);
 #endif
-    null_loglik = nrd.loglik;
+    //null_loglik = nrd.loglik;
 
     std::cout << " estimated null model ...";
 
@@ -191,7 +256,7 @@
 #if COXPH
     coxph_data rgd(phd,gtd,0,input_var.isIsInteractionExcluded());
 #else
-    regdata rgd(phd, gtd, 0,input_var.isIsInteractionExcluded());
+    regdata rgd(phd, gtd, 0, input_var.isIsInteractionExcluded());
 #endif
 
     std::cout << " formed regression object ...";
@@ -212,48 +277,7 @@
         {
             // open a file for output
             //_____________________
-
-            for (int i = 0; i < 5; i++)
-            {
-                outfile.push_back(new std::ofstream());
-            }
-
-            outfile[0]->open((outfilename_str + "_2df.out.txt").c_str());
-            outfile[1]->open((outfilename_str + "_add.out.txt").c_str());
-            outfile[2]->open((outfilename_str + "_domin.out.txt").c_str());
-            outfile[3]->open((outfilename_str + "_recess.out.txt").c_str());
-            outfile[4]->open((outfilename_str + "_over_domin.out.txt").c_str());
-
-            if (!outfile[0]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_2df.out.txt" << "\n";
-                exit(1);
-            }
-            if (!outfile[1]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_add.out.txt" << "\n";
-                exit(1);
-            }
-            if (!outfile[2]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_domin.out.txt" << "\n";
-                exit(1);
-            }
-            if (!outfile[3]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_recess.out.txt" << "\n";
-                exit(1);
-            }
-            if (!outfile[4]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_over_domin.out.txt" << "\n";
-                exit(1);
-            }
+            open_files_for_output(outfile, outfilename_str);
             //_____________________
 
             //Header
@@ -434,47 +458,7 @@
             //			outfilename_str="regression";
             //			}
 
-            for (int i = 0; i < 5; i++)
-            {
-                outfile.push_back(new std::ofstream());
-            }
-
-            outfile[0]->open((outfilename_str + "_2df.out.txt").c_str());
-            outfile[1]->open((outfilename_str + "_add.out.txt").c_str());
-            outfile[2]->open((outfilename_str + "_domin.out.txt").c_str());
-            outfile[3]->open((outfilename_str + "_recess.out.txt").c_str());
-            outfile[4]->open((outfilename_str + "_over_domin.out.txt").c_str());
-
-            if (!outfile[0]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_2df.out.txt" << "\n";
-                exit(1);
-            }
-            if (!outfile[1]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_add.out.txt" << "\n";
-                exit(1);
-            }
-            if (!outfile[2]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_domin.out.txt" << "\n";
-                exit(1);
-            }
-            if (!outfile[3]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_recess.out.txt" << "\n";
-                exit(1);
-            }
-            if (!outfile[4]->is_open())
-            {
-                std::cerr << "Cannot open file for writing: "
-                        << outfilename_str + "_over_domin.out.txt" << "\n";
-                exit(1);
-            }
+            open_files_for_output(outfile, outfilename_str);
         }
         else
         {
@@ -586,10 +570,11 @@
                     freq += snpdata1[ii] * 0.5;
                 }
         }
-        freq /= (double) gcount;
+        freq /= (double) (gcount);
         int poly = 1;
         if (fabs(freq) < 1.e-16 || fabs(1. - freq) < 1.e-16)
             poly = 0;
+
         if (fabs(mli.Rsq[csnp]) < 1.e-16)
             poly = 0;
 
@@ -758,7 +743,7 @@
                     //Oct 26, 2009
                     *chi2[model] << "nan";
                 }
-            } //end of moel cycle
+            } //end of model cycle
 
             //Han Chen
             *outfile[0] << beta_sebeta[0]->str() << input_var.getSep();
@@ -955,7 +940,6 @@
                 *outfile[0] << beta_sebeta[0]->str() << "\n";
             }
         }
-
         //clean chi2
         for (int i = 0; i < 5; i++)
         {
@@ -965,22 +949,8 @@
             //Oct 26, 2009
             chi2[i]->str("");
         }
+        update_progress_to_cmd_line(csnp, nsnps);
 
-        if (csnp % 1000 == 0)
-        {
-            if (csnp == 0)
-            {
-                fprintf(stdout, "Analysis: %6.2f ...",
-                        100. * double(csnp) / double(nsnps));
-            }
-            else
-            {
-                fprintf(stdout, "\b\b\b\b\b\b\b\b\b\b%6.2f ...",
-                        100. * double(csnp) / double(nsnps));
-            }
-            std::cout.flush();
-        }
-
     }
 
     fprintf(stdout, "\b\b\b\b\b\b\b\b\b\b%6.2f", 100.);



More information about the Genabel-commits mailing list