Implementation of MEDIAN() to calculate the median of a data sample in statistics.

Author: Joe7M

samples/distro-examples/tests/all.bas

@@ -217,6 +217,7 @@ print "SQR:" + SQR (x)
 print "SQUEEZE:" + SQUEEZE (s)
 print "STATMEAN:" + STATMEAN (1,2,3,4,5,6,7,8,9)
 print "STATMEANDEV:" + STATMEANDEV (1,2,3,4,5,6,7,8,9)
+print "STATMEDIAN:" + STATMEDIAN(1,2,3,4,5,6,7,8,9)
 print "STATSPREADP:" + STATSPREADP (1,2,3,4,5,6,7,8,9)
 print "STATSPREADS:" + STATSPREADS (1,2,3,4,5,6,7,8,9)
 print "STATSTD:" + STATSTD (1,2,3,4,5,6,7,8,9)

samples/distro-examples/tests/output/all.out

@@ -200,6 +200,7 @@ SQR:3.50713558335004
 SQUEEZE:catsanddogs
 STATMEAN:5
 STATMEANDEV:2.22222222222222
+STATMEDIAN:5
 STATSPREADP:6.66666666666667
 STATSPREADS:7.5
 STATSTD:2.73861278752583

src/common/blib_func.c

@@ -2620,6 +2620,7 @@ void cmd_genfunc(long funcCode, var_t *r) {
     //
     //
   case kwSTATMEANDEV:
+  case kwSTATMEDIAN:
   case kwSTATSTD:
   case kwSTATSPREADS:
   case kwSTATSPREADP:
@@ -2690,6 +2691,9 @@ void cmd_genfunc(long funcCode, var_t *r) {
       case kwSTATMEANDEV:
         r->v.n = statmeandev(dar, tcount);
         break;
+      case kwSTATMEDIAN:
+        r->v.n = statmedian(dar, tcount);
+        break;
       case kwSTATSTD:
         r->v.n = statstd(dar, tcount);
         break;

src/common/blib_math.c

@@ -232,6 +232,35 @@ var_num_t statmeandev(var_num_t *e, int count) {
   return sum / count;
 }
 
+//
+// compare function for median calculation
+//
+int median_compare (const void * a, const void * b) {
+  var_num_t fa = *(var_num_t*) a;
+  var_num_t fb = *(var_num_t*) b;
+  return (fa > fb) - (fa < fb);
+}
+
+//
+// Median
+//
+// Median is calculated using quicksort with a complexity of O(n log n).
+// Faster approch would be to implement i.e. quickselect with a 
+// complexity of O(n).
+//
+var_num_t statmedian(var_num_t *e, int count) {  
+  if (count == 0) {
+    return 0;
+  }
+
+  qsort(e, count, sizeof(var_num_t), median_compare);
+
+  if (count % 2 == 0)  
+      return (e[count/2] + e[count/2 - 1]) / 2;  
+  else  
+      return e[count/2];  
+}
+
 //
 // Standard deviation
 //

src/common/blib_math.h

@@ -75,6 +75,17 @@ var_num_t mat_determ(var_num_t *a, int n, double toler);
  */
 var_num_t statmeandev(var_num_t *e, int count);
 
+/**
+ * @ingroup math
+ * 
+ * Median
+ * 
+ * @param e array with numbers
+ * @param count number of elements of e
+ * @return the median
+ */
+var_num_t statmedian(var_num_t *e, int count);
+
 /**
  * @ingroup math
  * 

src/common/eval.c

@@ -1175,6 +1175,7 @@ static inline void eval_callf(var_t *r) {
   case kwSUMSV:
   case kwSTATMEAN:
   case kwSTATMEANDEV:
+  case kwSTATMEDIAN:
   case kwSTATSTD:
   case kwSTATSPREADS:
   case kwSTATSPREADP:

src/common/kw.h

@@ -348,6 +348,7 @@ enum func_keywords {
   kwSUM,
   kwSUMSV,
   kwSTATMEAN,
+  kwSTATMEDIAN,
   kwSTATSTD,
   kwSTATMEANDEV,
   kwSTATSPREADS,

src/languages/keywords.en.c

@@ -288,6 +288,7 @@ struct func_keyword_s func_table[] = {
 { "SUM",                        kwSUM },
 { "SUMSQ",                      kwSUMSV },
 { "STATMEAN",                   kwSTATMEAN },
+{ "STATMEDIAN",                 kwSTATMEDIAN},
 { "STATSTD",                    kwSTATSTD},
 { "STATMEANDEV",                kwSTATMEANDEV },
 { "STATSPREADS",                kwSTATSPREADS },