Changeset 502

Timestamp:

Jan 30, 2006, 11:05:20 PM (17 years ago)

Author:

Peter

Message:

improved median and percentile functions

Location:

trunk

Files:

• lib/classifier/ConsensusInputRanker.cc (modified) (1 diff)
• lib/statistics/WilcoxonFoldChange.cc (modified) (1 diff)
• lib/statistics/utility.h (modified) (3 diffs)
• test/statistics_test.cc (modified) (1 diff)

trunk/lib/classifier/ConsensusInputRanker.cc

@@ -43 +43 @@
         ranks[j]=input_rankers_[j].rank(i);
       }
-      medians[i].first = statistics::median_TMP(ranks);
+      medians[i].first = statistics::median(ranks);
       medians[i].second = i;
     }

trunk/lib/statistics/WilcoxonFoldChange.cc

@@ -35 +35 @@
     }
     if (absolute_)
-      return fabs(median_TMP(distance));
-    return median_TMP(distance);
+      return fabs(median(distance));
+    return median(distance);
   }
 

trunk/lib/statistics/utility.h

@@ -7 +7 @@
 
 #include <algorithm>
+#include <cassert>
 #include <cmath>
 #include <vector>
@@ -27 +28 @@
 
   ///
-  /// The percentile is determined by the \a p, a number between 0
-  /// and 100. The percentile is found by interpolation, using the
-  /// formula \f$ percentile = (1 - \delta) x_i + \delta x_{i+1} \f$
-  /// where \a p is floor\f$((n - 1)p/100)\f$ and \f$ \delta \f$ is
-  /// \f$ (n-1)p/100 - i \f$.Thus the minimum value of the vector is
-  /// given by p equal to zero, the maximum is given by p equal to
-  /// 100 and the median value is given by p equal to 50.
+  /// The percentile is determined by the \a p, a number between 0 and
+  /// 100. The percentile is found by interpolation, using the formula
+  /// \f$ percentile = (1 - \delta) x_i + \delta x_{i+1} \f$ where \a
+  /// p is floor\f$((n - 1)p/100)\f$ and \f$ \delta \f$ is \f$
+  /// (n-1)p/100 - i \f$.Thus the minimum value of the vector is given
+  /// by p equal to zero, the maximum is given by p equal to 100 and
+  /// the median value is given by p equal to 50. If @a sorted
+  /// is false (default), the vector is copied, the copy is sorted,
+  /// and then used to calculate the median.
   /// 
   /// @return \a p'th percentile
   ///
-  /// @note interface will change
-  ///
   template <class T>
-  double percentile_TMP(std::vector<T>& vec, double p)
+  double percentile(const std::vector<T>& vec, const double p, 
+                    const bool sorted=false)
   {
+    assert(!(p>100 && p<0));
+    if (sorted){
+      if (p>=100)
+        return vec.back();
+      double j = p/100 * (vec.size()-1);
+      int i = static_cast<int>(j);
+      return (1-j+floor(j))*vec[i] + (j-floor(j))*vec[i+1];
+    }
     if (p==100)
-      return  vec[vec.size()-1];
-    sort(vec.begin(), vec.end());
-    double j = p/100 * (vec.size()-1);
+      return  *std::max_element(vec.begin(),vec.end());
+    std::vector<T> v_copy(vec);
+    double j = p/100 * (v_copy.size()-1);
     int i = static_cast<int>(j);
-    return (1-j+floor(j))*vec[i] + (j-floor(j))*vec[i+1];
+    std::partial_sort(v_copy.begin(),v_copy.begin()+i+2 , v_copy.end());
+    return (1-j+floor(j))*v_copy[i] + (j-floor(j))*v_copy[i+1];
   
   }
 
   ///
-  /// Median is defined to be value in the middle. If number of
-  /// values is even median is the average of the two middle
-  /// values. @return median
+  /// Median is defined to be value in the middle. If number of values
+  /// is even median is the average of the two middle values.  the
+  /// median value is given by p equal to 50. If @a sorted is false
+  /// (default), the vector is copied, the copy is sorted, and then
+  /// used to calculate the median.
+  ///
+  /// @return median
   ///
   /// @note interface will change
   ///
   template <class T> 
-  inline double median_TMP(std::vector<T>& v) { return percentile_TMP(v,50.0); }
+  inline double median(const std::vector<T>& v, const bool sorted=false) 
+  { return percentile(v, 50.0, sorted); }
 
   ///
@@ -65 +81 @@
   /// is even median is the average of the two middle values. If @a
   /// sorted is true, the function assumes vector @a vec to be
-  /// sorted. If @a sorted is false, the vector is copied, the copy
-  /// is sorted, and then used to calculate the median.
+  /// sorted. If @a sorted is false, the vector is copied, the copy is
+  /// sorted (default), and then used to calculate the median.
   ///
   /// @return median
   ///
-  double median(const gslapi::vector& vec, const bool sorted=true);
+  double median(const gslapi::vector& vec, const bool sorted=false);
 
 }} // of namespace statistics and namespace theplu

trunk/test/statistics_test.cc

@@ -15 +15 @@
   for (unsigned int i=0; i<10; i++)
     data.push_back(static_cast<double>(i));
-  double m=theplu::statistics::median_TMP(data);
+  double m=theplu::statistics::median(data);
   if (m!=4.5)
     return -1;