Changeset 1062

Timestamp:

May 13, 2009, 8:09:05 AM (14 years ago)

Author:

Jari Häkkinen

Message:

Addresses #118. Added support to use either mean or median in target calculation.

File:

• plugins/base2/net.sf.basedb.normalizers/trunk/src/c++/bin/qQN.cc (modified) (6 diffs)

plugins/base2/net.sf.basedb.normalizers/trunk/src/c++/bin/qQN.cc

@@ -35 +35 @@
 #include <yat/utility/OptionSwitch.h>
 #include <yat/utility/stl_utility.h>
+#include <yat/statistics/utility.h>
 
 #include <cstdlib>
@@ -45 +46 @@
 
 
-void create_target(std::vector<double>&, const MatrixWeighted&);
+void create_target(std::vector<double>&, const MatrixWeighted&, bool use_median);
 void create_target(std::vector<double>&, const MatrixWeighted&,
-                   const std::string&);
+                   const std::string&, bool use_median);
 
 /**
@@ -62 +63 @@
   OptionInFile indata(cmd, "in-data", "data to be normalized");
   OptionOutFile outdata(cmd, "out-data", "normalized data");
+  OptionSwitch target_median(cmd, "median",
+                             "use median instead of mean for targets", false);
   OptionHelp help(cmd);
   help.synopsis()=(std::string("See ") +
@@ -103 +106 @@
 
   std::vector<double> target;
-  ( assay.present() ? create_target(target,m,assay.value()) :
-                      create_target(target,m) );
+  ( assay.present() ?
+    create_target(target, m, assay.value(), target_median.value()) :
+    create_target(target, m, target_median.value()) );
   std::transform(target.begin(), target.end(),
                  target.begin(), theplu::yat::utility::Log<double>());
@@ -135 +139 @@
 
 void create_target(std::vector<double>& t, const MatrixWeighted& m,
-                   const std::string& assay)
-{
+                   const std::string& assay, bool use_median)
+{
+  std::vector< std::vector<double> > calc_support(m.rows());
+  for (std::vector< std::vector<double> >::iterator i=calc_support.begin();
+       i!=calc_support.end(); ++i)
+    i->reserve(m.columns());
+
   std::ifstream is(assay.c_str());
   std::string line;
   size_t column=0;
-  std::vector<size_t> column_contribs(m.rows(),0);
-  std::vector<double> temp_target(m.rows(),0.0);
   while (getline(is, line)) {
     size_t found=line.find("yes");
-    if (found!=std::string::npos) {
-      // use this assay as a part of reference
+    if (found!=std::string::npos) // use this assay as a part of reference
       for (size_t row=0; row<m.rows(); ++row)
-        if (m(row,column).weight()) { // weight either 0 or 1
-          temp_target[row]+=m(row,column).data();
-          ++column_contribs[row];
-        }
-    }
+        if (m(row,column).weight())       // weight either 0 or 1
+          calc_support[row].push_back(m(row,column).data());
     ++column;
     if (column>m.columns())
       throw std::runtime_error("Too many annotation columns wrt data matrix");
   }
+
   t.reserve(m.rows());
-  for (size_t row=0; row<m.rows(); ++row)
-    if (column_contribs[row])
-      t.push_back(temp_target[row]/=column_contribs[row]);
+  for (size_t row=0; row<m.rows(); ++row) {
+    const std::vector<double>& cs=calc_support[row];
+    if (!cs.empty())
+      t.push_back( use_median ?
+                   theplu::yat::statistics::median(cs.begin(), cs.end()) :
+                   std::accumulate(cs.begin(), cs.end(), 0.0) / cs.size() );
+  }
   if (!t.size())
     throw std::runtime_error("Not a well defined reference, aborting");
@@ -165 +173 @@
 
 
-void create_target(std::vector<double>& t, const MatrixWeighted& m)
-{
-  std::vector<double> temp_target(m.rows(),0.0);
+void create_target(std::vector<double>& t, const MatrixWeighted& m,
+                   bool use_median)
+{
   t.reserve(m.rows());
+  std::vector<double> cs;
+  cs.reserve(m.columns());
   for (size_t row=0; row<m.rows(); ++row) {
-    size_t column_contribs=0;
+    cs.clear();
     for (size_t column=0; column<m.columns(); ++column)
-      if (m(row,column).weight()) { // weight either 0 or 1
-        temp_target[row]+=m(row,column).data();
-        ++column_contribs;
-      }
-    if (column_contribs)
-      t.push_back(temp_target[row]/=column_contribs);
+      if (m(row,column).weight())      // weight either 0 or 1
+        cs.push_back(m(row,column).data());
+    if (!cs.empty())
+      t.push_back( use_median ?
+                   theplu::yat::statistics::median(cs.begin(), cs.end()) :
+                   std::accumulate(cs.begin(), cs.end(), 0.0) / cs.size() );
   }
   if (!t.size())