[Genabel-commits] r1820 - in pkg/OmicABELnoMM: examples examples/small src tests

noreply at r-forge.r-project.org noreply at r-forge.r-project.org
Tue Sep 9 17:43:18 CEST 2014


Author: afrank
Date: 2014-09-09 17:43:18 +0200 (Tue, 09 Sep 2014)
New Revision: 1820

Added:
   pkg/OmicABELnoMM/examples/small/
   pkg/OmicABELnoMM/examples/small/XL.fvd
   pkg/OmicABELnoMM/examples/small/XL.fvi
   pkg/OmicABELnoMM/examples/small/XR.fvd
   pkg/OmicABELnoMM/examples/small/XR.fvi
   pkg/OmicABELnoMM/examples/small/Y.fvd
   pkg/OmicABELnoMM/examples/small/Y.fvi
   pkg/OmicABELnoMM/examples/small/bpre.fvd
Modified:
   pkg/OmicABELnoMM/examples/dosages_2.txt
   pkg/OmicABELnoMM/src/AIOwrapper.cpp
   pkg/OmicABELnoMM/src/AIOwrapper.h
   pkg/OmicABELnoMM/src/main.cpp
   pkg/OmicABELnoMM/tests/test.cpp
Log:
Fixed Modeling Bug. Changed Help Cmd. Improved performance for Std. Models where factors are zero. Added smaller file examples.

Modified: pkg/OmicABELnoMM/examples/dosages_2.txt
===================================================================
--- pkg/OmicABELnoMM/examples/dosages_2.txt	2014-09-09 13:54:05 UTC (rev 1819)
+++ pkg/OmicABELnoMM/examples/dosages_2.txt	2014-09-09 15:43:18 UTC (rev 1820)
@@ -1 +1 @@
-2 1
\ No newline at end of file
+1 1
\ No newline at end of file

Added: pkg/OmicABELnoMM/examples/small/XL.fvd
===================================================================
(Binary files differ)


Property changes on: pkg/OmicABELnoMM/examples/small/XL.fvd
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/OmicABELnoMM/examples/small/XL.fvi
===================================================================
(Binary files differ)


Property changes on: pkg/OmicABELnoMM/examples/small/XL.fvi
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/OmicABELnoMM/examples/small/XR.fvd
===================================================================
(Binary files differ)


Property changes on: pkg/OmicABELnoMM/examples/small/XR.fvd
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/OmicABELnoMM/examples/small/XR.fvi
===================================================================
(Binary files differ)


Property changes on: pkg/OmicABELnoMM/examples/small/XR.fvi
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/OmicABELnoMM/examples/small/Y.fvd
===================================================================
(Binary files differ)


Property changes on: pkg/OmicABELnoMM/examples/small/Y.fvd
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/OmicABELnoMM/examples/small/Y.fvi
===================================================================
(Binary files differ)


Property changes on: pkg/OmicABELnoMM/examples/small/Y.fvi
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Added: pkg/OmicABELnoMM/examples/small/bpre.fvd
===================================================================
(Binary files differ)


Property changes on: pkg/OmicABELnoMM/examples/small/bpre.fvd
___________________________________________________________________
Added: svn:executable
   + *
Added: svn:mime-type
   + application/octet-stream

Modified: pkg/OmicABELnoMM/src/AIOwrapper.cpp
===================================================================
--- pkg/OmicABELnoMM/src/AIOwrapper.cpp	2014-09-09 13:54:05 UTC (rev 1819)
+++ pkg/OmicABELnoMM/src/AIOwrapper.cpp	2014-09-09 15:43:18 UTC (rev 1820)
@@ -32,7 +32,7 @@
 
 
     Fhandler->use_dosages = params.dosages;
-    if(params.dosages && Fhandler->model ==-1)
+    if(params.dosages && params.model ==-1)
     {
         cout << "Requested dosages model wihtout a valid model!" << endl;
         exit(1);
@@ -40,6 +40,9 @@
     Fhandler->not_done = true;
     Fhandler->model = params.model;
     Fhandler->fname_dosages = params.fname_dosages;
+    Fhandler->dosage_skip = 0;
+    Fhandler->fileR = 0;
+    Fhandler->add_dosages = false;
 
 
 
@@ -91,26 +94,10 @@
 
 
 
-        int Aname_idx=params.n*ARfvi->fvi_header.namelength;//skip the names of the rows
-        if(Fhandler->use_dosages)
-        {
-            for(int i = 0; i < params.m; i++)
-            {
-                Fhandler->ARnames.push_back(string(&(ARfvi->fvi_data[Aname_idx])));
-                Aname_idx += ARfvi->fvi_header.namelength*Fhandler->fileR;
-            }
-        }
-        else
-        {
-            for(int i = 0; i < params.m*params.r; i++)
-            {
-                Fhandler->ARnames.push_back(string(&(ARfvi->fvi_data[Aname_idx])));
-                Aname_idx += ARfvi->fvi_header.namelength;
-            }
-        }
 
 
-        Aname_idx=params.n*ALfvi->fvi_header.namelength;
+
+        int Aname_idx=params.n*ALfvi->fvi_header.namelength;
         for(int i = 0; i < params.l; i++)
         {
             Fhandler->ALnames.push_back(string(&(ALfvi->fvi_data[Aname_idx])));
@@ -174,37 +161,47 @@
 
         break;
         case 0://add
-            if(Fhandler->fileR != 3)
+            if(Fhandler->fileR != 3 && Fhandler->fileR != 2)
             {
-                cout << "The amount of columns per Independent Variable (--ngpred) is not 3 for a valid Additive Model!" << endl;
+                cout << "The amount of columns per Independent Variable (--ngpred) is not 2 or 3 for a valid Additive Model!" << endl;
                 exit(1);
             }
-            Fhandler->dosages[0] = 2;Fhandler->dosages[1] = 1;Fhandler->dosages[2] = 0;
+            Fhandler->dosages[0] = 2;Fhandler->dosages[1] = 1;//Fhandler->dosages[2] = 0;
+            if(Fhandler->fileR == 3)
+                Fhandler->dosage_skip = 1;//last columns is irrelevant
             params.r = 1;
             Fhandler->add_dosages = true;
         break;
         case 1://dom
-            if(Fhandler->fileR != 3)
+            if(Fhandler->fileR != 3 && Fhandler->fileR != 2)
             {
-                cout << "The amount of columns per Independent Variable (--ngpred) is not 3 for a valid Dominant Model!" << endl;
+                cout << "The amount of columns per Independent Variable (--ngpred) is not 2 or 3 for a valid Dominant Model!" << endl;
                 exit(1);
             }
-            Fhandler->dosages[0] = 1;Fhandler->dosages[1] = 1;Fhandler->dosages[2] = 0;
+            Fhandler->dosages[0] = 1;Fhandler->dosages[1] = 1;//Fhandler->dosages[2] = 0;
+            if(Fhandler->fileR == 3)
+                Fhandler->dosage_skip = 1;//last columns is irrelevant
             params.r = 1;
             Fhandler->add_dosages = true;
         break;
         case 2://rec
-            if(Fhandler->fileR != 3)
+            if(Fhandler->fileR != 3 && Fhandler->fileR != 2)
             {
-                cout << "The amount of columns per Independent Variable (--ngpred) is not 3 for a valid Recessive Model!" << endl;
+                cout << "The amount of columns per Independent Variable (--ngpred) is not 2 or 3 for a valid Recessive Model!" << endl;
                 exit(1);
             }
-            Fhandler->dosages[0] = 1;Fhandler->dosages[1] = 0;Fhandler->dosages[2] = 0;
+            Fhandler->dosages[0] = 1;//Fhandler->dosages[1] = 0;Fhandler->dosages[2] = 0;
+            Fhandler->dosage_skip = 2;
+            if(Fhandler->fileR == 3)
+                (Fhandler->dosage_skip)--;//last 2 columns is irrelevant
             params.r = 1;
             Fhandler->add_dosages = true;
         break;
         case 3://linear
+
             read_dosages(params.fname_dosages,Fhandler->fileR,Fhandler->dosages);
+            Fhandler->add_dosages = false;
+            Fhandler->dosage_skip = 0 ;
         break;
         case 4://additive
             read_dosages(params.fname_dosages,Fhandler->fileR,Fhandler->dosages);
@@ -212,13 +209,35 @@
             Fhandler->add_dosages = true;
         break;
         }
-    }
+    }
+
+
+
 
+
     params.p = params.l + params.r;
 
     if(!Fhandler->fakefiles)
     {
 
+        int Aname_idx=params.n*ARfvi->fvi_header.namelength;//skip the names of the rows
+        if(Fhandler->use_dosages && Fhandler->add_dosages)
+        {
+            for(int i = 0; i < params.m; i++)
+            {
+                Fhandler->ARnames.push_back(string(&(ARfvi->fvi_data[Aname_idx])));
+                Aname_idx += ARfvi->fvi_header.namelength*Fhandler->fileR;
+            }
+        }
+        else
+        {
+            for(int i = 0; i < params.m*params.r; i++)
+            {
+                Fhandler->ARnames.push_back(string(&(ARfvi->fvi_data[Aname_idx])));
+                Aname_idx += ARfvi->fvi_header.namelength;
+            }
+        }
+
         if(params.l > 255 || params.p > 255 || params.n > 65535)//can remove if fixed for output files
         {
             cout << "Warning, output binary format does not yet support current problem sizes for the provided p, l, r and n." << endl;
@@ -245,11 +264,11 @@
         fp_InfoResults.write( (char*)&params.storePInd,sizeof(bool));
 
 
-        int Aname_idx=(params.n+1)*ALfvi->fvi_header.namelength;//skip the names of the rows + intercept
+        Aname_idx=(params.n+1)*ALfvi->fvi_header.namelength;//skip the names of the rows + intercept
         fp_InfoResults.write( (char*)&ALfvi->fvi_data[Aname_idx],ALfvi->fvi_header.namelength*(params.l-1)*sizeof(char));
 
         Aname_idx=params.n*ARfvi->fvi_header.namelength;//skip the names of the rows
-        if(Fhandler->use_dosages)
+        if(Fhandler->use_dosages && Fhandler->add_dosages)
         {
             for(int i = 0; i < params.m; i++)
             {
@@ -286,9 +305,8 @@
 
 
 
+    //
 
-
-
 }
 
 void AIOwrapper::finalize()
@@ -540,52 +558,53 @@
             }
             else
             {
+                //cout << " "<< Fhandler->use_dosages << " " <<  Fhandler->add_dosages <<" " <<  Fhandler->model << " " << Fhandler->fileR << " " << Fhandler->dosage_skip << endl;
 
+
                 list< pair<int,int> >* excl_List = Fhandler->excl_List;
 
                 int chunk_size_buff;
-                int buff_pos=0;
+
                 int file_pos;
-                float* destination = Fhandler->ArDosage;
+                int buff_pos = 0;
 
                 if(Fhandler->use_dosages)
                 {
 
-                    if(!Fhandler->add_dosages)
-                    {
-                        destination = tobeFilled->buff;//no need to use temp variable
-                    }
 
                     for(int i = 0; i < tmp_ar_blockSize; i++)
                     {
                         buff_pos=0;
-                        for(int ii = 0; ii < Fhandler->fileR; ii++)
+                        for(int ii = 0; ii < (Fhandler->fileR-Fhandler->dosage_skip); ii++)
                         {
+
                             for (list<  pair<int,int>  >::iterator it=excl_List->begin(); it != excl_List->end(); ++it)
                             {
-                                file_pos = Fhandler->fileN*i + it->first;
+                                file_pos = i*(Fhandler->fileN*Fhandler->fileR)+ ii*Fhandler->fileN + it->first;
                                 fseek ( fp_Ar , file_pos*sizeof(type_precision) , SEEK_SET );
 
                                 chunk_size_buff = it->second;
 
-                                size_t result = fread (&(destination[buff_pos]),sizeof(type_precision),chunk_size_buff,fp_Ar); result++;
+                                size_t result = fread (&(Fhandler->ArDosage[buff_pos]),sizeof(type_precision),chunk_size_buff,fp_Ar); result++;
                                 buff_pos += chunk_size_buff;
                             }
                         }
 
                         if(Fhandler->add_dosages)
                         {
+                            //cout << "adding ";
                             cblas_sgemm(CblasColMajor, CblasNoTrans, CblasNoTrans,
-                                Fhandler->n, 1, Fhandler->fileR, 1.0, Fhandler->ArDosage, Fhandler->n, Fhandler->dosages,Fhandler->fileR ,
+                                Fhandler->n, 1, Fhandler->fileR-Fhandler->dosage_skip, 1.0, Fhandler->ArDosage, Fhandler->n, Fhandler->dosages,Fhandler->fileR-Fhandler->dosage_skip ,
                                     0.0, &(tobeFilled->buff[i*Fhandler->n]), Fhandler->n);
                         }
                         else
                         {
+                            //cout << "mult ";
                             for(int ii = 0; ii < Fhandler->fileR; ii++)
                             {
                                 for(int k=0; k < Fhandler->n; k++)
                                 {
-                                    destination[Fhandler->n*ii+k] *= Fhandler->dosages[ii];
+                                   tobeFilled->buff[i*Fhandler->n*Fhandler->r + ii*Fhandler->n+k] = Fhandler->ArDosage[Fhandler->n*ii+k] * Fhandler->dosages[ii];
                                 }
                             }
                         }
@@ -608,7 +627,6 @@
                             size_t result = fread (&tobeFilled->buff[buff_pos],sizeof(type_precision),chunk_size_buff,fp_Ar); result++;
                             buff_pos += chunk_size_buff;
 
-
                         }
                     }
                 }

Modified: pkg/OmicABELnoMM/src/AIOwrapper.h
===================================================================
--- pkg/OmicABELnoMM/src/AIOwrapper.h	2014-09-09 13:54:05 UTC (rev 1819)
+++ pkg/OmicABELnoMM/src/AIOwrapper.h	2014-09-09 15:43:18 UTC (rev 1820)
@@ -76,8 +76,9 @@
     int l;
     int p;
 
-    int model;
+    int model;
 
+
     int Ar_Amount;
     int Ar_blockSize;
     int Ar_to_readSize;
@@ -92,7 +93,8 @@
     bool not_done;
     bool reset_wait;
     bool use_dosages;
-    bool add_dosages;
+    bool add_dosages;
+    int dosage_skip;
 
     int seed;
     int Aseed;

Modified: pkg/OmicABELnoMM/src/main.cpp
===================================================================
--- pkg/OmicABELnoMM/src/main.cpp	2014-09-09 13:54:05 UTC (rev 1819)
+++ pkg/OmicABELnoMM/src/main.cpp	2014-09-09 15:43:18 UTC (rev 1820)
@@ -236,35 +236,35 @@
             break;
 
         case 'j':
-            params.model = 1;
+            params.model = 0;
             params.dosages = true;
 
             cout << "-j Using Additive Model with (2*AA,1*AB,0*BB) effects"<< endl;
             break;
 
         case 'k':
-            params.model = 2;
+            params.model = 1;
             params.dosages = true;
 
             cout << "-j Using Dominant Model with (1*AA,1*AB,0*BB) effects"<< endl;
             break;
 
         case 'l':
-            params.model = 3;
+            params.model = 2;
             params.dosages = true;
 
-            cout << "-j Using Recessive Model with (0*AA,0*AB,1*BB) effects"<< endl;
+            cout << "-j Using Recessive Model with (1*AA,0*AB,0*BB) effects"<< endl;
             break;
 
         case 'z':
-            params.model = 4;
+            params.model = 3;
             params.dosages = true;
 
             cout << "-z Using Custom Linear Model with parameters read from the file "<< params.fname_dosages << endl;
             break;
 
         case 'y':
-            params.model = 5;
+            params.model = 4;
             params.dosages = true;
 
             cout << "-z Using Custom Additive Model with parameters read from the file "<< params.fname_dosages << endl;

Modified: pkg/OmicABELnoMM/tests/test.cpp
===================================================================
--- pkg/OmicABELnoMM/tests/test.cpp	2014-09-09 13:54:05 UTC (rev 1819)
+++ pkg/OmicABELnoMM/tests/test.cpp	2014-09-09 15:43:18 UTC (rev 1820)
@@ -117,9 +117,9 @@
     params.fnameAR="examples/XR";
     params.fnameY="examples/Y";
     params.fnameOutFiles="resultsSig";
-//    params.dosages = true;
-//    params.model = 4;
-//    params.fname_dosages = "examples/dosages_2.txt";
+    params.dosages = true;
+    params.model = 0;
+    params.fname_dosages = "examples/dosages_2.txt";
 
 
     for(int th = 0; th < max_threads; th++)



More information about the Genabel-commits mailing list