Changeset 1536

Timestamp:

Sep 25, 2008, 9:31:57 PM (15 years ago)

Author:

Peter

Message:

adding a Gauss normalizer

Location:

trunk

Files:

• test/normalization_test.cc (modified) (4 diffs)
• yat/normalizer/Gauss.h (copied) (copied from trunk/yat/normalizer/Spearman.h) (5 diffs)
• yat/normalizer/Makefile.am (modified) (1 diff)

trunk/test/normalization_test.cc

@@ -24 +24 @@
 #include "yat/normalizer/Centralizer.h"
 #include "yat/normalizer/ColumnNormalizer.h"
-#include "yat/utility/DataIterator.h"
+#include "yat/normalizer/Gauss.h"
 #include "yat/normalizer/QuantileNormalizer.h"
 #include "yat/normalizer/RowNormalizer.h"
 #include "yat/normalizer/Spearman.h"
 #include "yat/normalizer/Zscore.h"
+#include "yat/utility/DataIterator.h"
 
 #include "yat/utility/DataIterator.h"
@@ -41 +42 @@
 void test_centralizer(test::Suite&);
 void test_column_normalize(test::Suite&);
+void test_gauss_normalize(test::Suite&);
 void test_quantile_normalize(test::Suite&);
 void test_row_normalize(test::Suite&);
@@ -55 +57 @@
   test_column_normalize(suite);
   test_quantile_normalize(suite);
+  test_gauss_normalize(suite);
   test_row_normalize(suite);
   test_spearman(suite);
@@ -178 +181 @@
 
 
+void test_gauss_normalize(test::Suite& suite)
+{
+  suite.err() << "Testing Gauss\n";
+  normalizer::Gauss gauss;
+  std::vector<double> vec;
+  vec.push_back(1);
+  gauss(vec.begin(), vec.end(), vec.begin());
+  suite.add(suite.equal(vec.front(), 0));
+  vec.push_back(1);
+  gauss(vec.begin(), vec.end(), vec.begin());
+  suite.add(suite.equal(vec.front(), -vec.back()));
+
+}
+
 void test_spearman_weighted(test::Suite& suite)
 {

trunk/yat/normalizer/Gauss.h

@@ -1 @@
-#ifndef _theplu_yat_normalizer_spearman_
-#define _theplu_yat_normalizer_spearman_
+#ifndef _theplu_yat_normalizer_gauss_
+#define _theplu_yat_normalizer_gauss_
 
 /*
@@ -21 +21 @@
 */
 
-#include "yat/utility/DataIterator.h"
+#include "Spearman.h"
 #include "yat/utility/iterator_traits.h"
-#include "yat/utility/sort_index.h"
-#include "yat/utility/WeightIterator.h"
 
-#include <algorithm>
-#include <functional>
-#include <vector>
+#include <gsl/gsl_cdf.h>
 
 namespace theplu {
@@ -35 +31 @@
 
   /**
-     \brief Replace elements with normalized rank
+     \brief Gaussian Normalizer
 
      \since New in yat 0.5
    */
-  class Spearman
+  class Gauss
   {
   public:
-    /**
-       \brief default constructor
-     */
-    Spearman(void){}
-
     /**
        It is possible to centralize a range "in place"; it is
@@ -52 +43 @@
        same.
 
-       Each element x is replaced by \f$ \frac{\sum I(x_i-x)
-       w_i}{\sum w_i} \f$ where I(x) = 1 for x>0, I(x) = 0.5 for x=0,
-       and I(x) = 0 for x<0.
+       The range is first rank normalized using Spearman, after which
+       each element is between 0 and unity. Second each element is
+       replaced by inverse cumulative standard Gaussian distribution.
+       
+       After normalization the range will follow a standard Gaussian
+       distribution (mean zero and unity variance).
+
+       \see gsl_cdf_ugaussian_Pinv
 
        \return result + (last-first)
@@ -62 +58 @@
                                     RandomAccessIterator result) const
     {
-      typename utility::weighted_iterator_traits<ForwardIterator>::type tag;
-      return normalize(first, last, result, tag);
+      Spearman spearman;
+      spearman(first, last, result);
+      RandomAccessIterator end = result + (last-first);
+      utility::iterator_traits<RandomAccessIterator> trait;
+      while (result != end) {
+        trait.data(result) = gsl_cdf_ugaussian_Pinv(trait.data(result));
+        ++result;
+      }
+      
+      return result;
     }
 
 
   private:
-    // unweighted version
-    template<typename ForwardIterator, typename RandomAccessIterator>
-    RandomAccessIterator normalize(ForwardIterator first, ForwardIterator last,
-                                   RandomAccessIterator result, 
-                                   utility::unweighted_iterator_tag) const
-    {
-      std::vector<size_t> perm;
-      utility::sort_index(first, last, perm);
-      double n = perm.size();
-      size_t i=0; 
-      while ( i<perm.size() ) {
-        size_t min_i = i;
-        while (i<perm.size() && first[perm[i]]<=first[perm[min_i]])
-          ++i;
-        double res = static_cast<double>(i + min_i)/(2*n);
-        for ( ; min_i < i; ++min_i)
-          result[perm[min_i]] = res;
-      }
-      return result + std::distance(first, last);
-    }
-
-
-    // weighted version
-    template<typename ForwardIterator, typename RandomAccessIterator>
-    RandomAccessIterator normalize(ForwardIterator first, ForwardIterator last,
-                                   RandomAccessIterator result, 
-                                   utility::weighted_iterator_tag) const
-    {
-      std::copy(utility::weight_iterator(first), 
-                utility::weight_iterator(last),
-                utility::weight_iterator(result));
-      // set values with w=0 to 0 to avoid problems with NaNs
-      utility::iterator_traits<ForwardIterator> trait;
-      for (ForwardIterator i=first; i!=last; ++i)
-        if (trait.weight(i)==0)
-          trait.data(i)=0.0;
-
-      std::vector<size_t> perm(std::distance(first, last));
-      utility::sort_index(utility::data_iterator(first), 
-                          utility::data_iterator(last), perm);
-      utility::iterator_traits<RandomAccessIterator> rtrait;
-
-      double sum_w=0;
-      size_t i=0; 
-      while ( i<perm.size() ) {
-        double w=0;
-        size_t min_i = i;
-        while (i<perm.size() && (trait.weight(first+perm[i])==0 ||
-                                 trait.data(first+perm[i]) <= 
-                                 trait.data(first+perm[min_i])) ) {
-          w += trait.weight(first+perm[i]);
-          ++i;
-        }
-        double res=sum_w + 0.5*w;
-        for ( size_t j=min_i; j<i; ++j)
-          rtrait.data(result+perm[j]) = res;
-        sum_w += w;
-      }       
-
-      size_t n = std::distance(first, last);
-      std::transform(utility::data_iterator(result),
-                     utility::data_iterator(result+n),
-                     utility::data_iterator(result),
-                     std::bind2nd(std::divides<double>(), sum_w));
-      return result + n;
-    }
-
   };
 

trunk/yat/normalizer/Makefile.am

@@ -26 +26 @@
 
 include_normalizer_HEADERS = Centralizer.h ColumnNormalizer.h \
-  QuantileNormalizer.h RowNormalizer.h Spearman.h Zscore.h
+  Gauss.h QuantileNormalizer.h RowNormalizer.h Spearman.h Zscore.h