home *** CD-ROM | disk | FTP | other *** search
/ io Programmo 41 / IOPROG_41.ISO / soft / c++ / NUMCPP11.ZIP / descript.hpp < prev    next >
Encoding:
C/C++ Source or Header  |  1999-05-15  |  12.5 KB  |  397 lines

  1. //===================================================================
  2. // descript.hpp
  3. //
  4. // Version 1.1
  5. //
  6. // Written by:
  7. //   Brent Worden
  8. //   WordenWare
  9. //   email:  Brent@Worden.org
  10. //
  11. // Copyright (c) 1998-1999 WordenWare
  12. //
  13. // Created:  August 28, 1998
  14. // Revised:  April 10, 1999
  15. //===================================================================
  16.  
  17. #ifndef _DESCRIPT_HPP_
  18. #define _DESCRIPT_HPP_
  19.  
  20. #include <cmath>
  21.  
  22. #include <algorithm>
  23. #include <numeric>
  24.  
  25. #include "algorthm.hpp"
  26. #include "residual.hpp"
  27. #include "vector.hpp"
  28.  
  29. NUM_BEGIN
  30.  
  31. template<class Iter>
  32. double centmom(Iter first, Iter last, double n)
  33. //-------------------------------------------------------------------
  34. // Returns the n-th central moment of the elements in [first, last).
  35. //-------------------------------------------------------------------
  36. {
  37.     Vector<double> r(length(first, last));
  38.     std::copy(first, last, r.begin());
  39.     
  40.     nresid(r.begin(), r.end(), mean(first, last), n);
  41.     return mean(r.begin(), r.end());
  42. };
  43.  
  44. template<class Iter>
  45. inline double coeffvar(Iter first, Iter last)
  46. //-------------------------------------------------------------------
  47. // Returns the coefficient of variation of the elements in
  48. // [first, last).
  49. //-------------------------------------------------------------------
  50. {
  51.     return stddev(first, last) / mean(first, last) * 100.0;
  52. };
  53.  
  54. template<class Iter>
  55. double geomean(Iter first, Iter last)
  56. //-------------------------------------------------------------------
  57. // Returns the geometric mean of the elements in [first, last).
  58. //-------------------------------------------------------------------
  59. {
  60.     double p = 1.0 / length(first, last);
  61.     double m = 1.0;
  62.     
  63.     while(first != last){
  64.         m *= pow(*first, p);
  65.         ++first;
  66.     };
  67.     
  68.     return m;
  69. };
  70.  
  71. template<class Iter>
  72. double harmean(Iter first, Iter last)
  73. //-------------------------------------------------------------------
  74. // Returns the harmonic mean of the elements in [first, last).
  75. //-------------------------------------------------------------------
  76. {
  77.     double n = length(first, last);
  78.     double m = 0.0;
  79.     
  80.     while(first != last){
  81.         if(*first){
  82.             m += 1.0 / *first;
  83.         }
  84.         ++first;
  85.     }
  86.     
  87.     return m;
  88. };
  89.  
  90. template<class Iter>
  91. inline double iqr(Iter first, Iter last)
  92. //-------------------------------------------------------------------
  93. // Returns the interquartile range of the elements in [first, last).
  94. //-------------------------------------------------------------------
  95. {
  96.     return quartile3(first, last) - quartile1(first, last);
  97. };
  98.  
  99. template<class Iter>
  100. inline double kurtosis(Iter first, Iter last)
  101. //-------------------------------------------------------------------
  102. // Returns the kurtosis of the elements in [first, last).
  103. //-------------------------------------------------------------------
  104. {
  105.     return centmom(first, last, 4) / pow(variance(first, last), 2);
  106. };
  107.  
  108. template<class Iter>
  109. int length(Iter first, Iter last)
  110. //-------------------------------------------------------------------
  111. // Returns the number of elements in [first, last).
  112. //-------------------------------------------------------------------
  113. {
  114.     int len = 0;
  115.     while(first != last){
  116.         ++first;
  117.         ++len;
  118.     }
  119.     return len;
  120. };
  121.  
  122. template<class Iter>
  123. double mean(Iter first, Iter last)
  124. //-------------------------------------------------------------------
  125. // Returns the mean of the elements in [first, last).
  126. //-------------------------------------------------------------------
  127. {
  128.     double m = 0.0;
  129.     double n = 0.0;
  130.     while(first != last){
  131.         n += 1.0;
  132.         m += (*first - m) / n;
  133.         ++first;
  134.     }
  135.     return m;
  136. };
  137.  
  138. template<class Iter>
  139. double meandev(Iter first, Iter last)
  140. //-------------------------------------------------------------------
  141. // Returns the mean deviation of the elements in [first, last).
  142. //-------------------------------------------------------------------
  143. {
  144.     Vector<double> r(length(first, last));
  145.     std::copy(first, last, r.begin());
  146.     aresid(r.begin(), r.end(), mean(first, last));
  147.     return mean(r.begin(), r.end());
  148. };
  149.  
  150. template<class Iter>
  151. double meddev(Iter first, Iter last)
  152. //-------------------------------------------------------------------
  153. // Returns the median deviation of the elements in [first, last).
  154. //-------------------------------------------------------------------
  155. {
  156.     Vector<double> r(length(first, last));
  157.     std::copy(first, last, r.begin());
  158.     aresid(r.begin(), r.end(), median(first, last));
  159.     return median(r.begin(), r.end());
  160. };
  161.  
  162. template<class Iter>
  163. double median(Iter first, Iter last)
  164. //-------------------------------------------------------------------
  165. // Returns the median of the elements in [first, last).
  166. //-------------------------------------------------------------------
  167. {
  168.     int n = length(first, last);
  169.     Vector<double> v(n);
  170.     std::copy(first, last, v.begin());
  171.     Vector<double>::iterator nth = v.begin() + n / 2;
  172.     std::nth_element(v.begin(), nth, v.end());
  173.     
  174.     if(n % 2 == 0){ // even n
  175.         double b = *nth;
  176.         nth = v.begin() + n / 2 - 1;
  177.         std::nth_element(v.begin(), nth, v.end());
  178.         return average(*nth, b);
  179.     } else {
  180.         return *nth;
  181.     }
  182. };
  183.  
  184. template<class Iter>
  185. double quantile(Iter first, Iter last, double q)
  186. //-------------------------------------------------------------------
  187. // Returns the q*100 percentile of the elements in [first, last).
  188. //-------------------------------------------------------------------
  189. {
  190.     int n = length(first, last);
  191.     Vector<double> vec(n);
  192.     std::copy(first, last, vec.begin());
  193.     double u, f, v, ans;
  194.     int i;
  195.     
  196.     if(q < 0.0){
  197.         q = 0.0;
  198.     }
  199.     if(q > 1.0){
  200.         q = 1.0;
  201.     }
  202.     
  203.     std::sort(vec.begin(), vec.end());
  204.     if(q < (v = 1.0 / (2.0 * n))){
  205.         ans = *(vec.begin());
  206.     } else if(q > 1.0 - v){
  207.         ans = *(vec.end() - 1);
  208.     } else {
  209.         u = n * q + .5;
  210.         i = int(u);
  211.         f = u - i;
  212.         ans = (1.0 - f) * *(vec.begin() + i - 1) + f * *(vec.begin() + i);
  213.     }
  214.     
  215.     return ans;
  216. };
  217.  
  218. template<class Iter>
  219. inline double quartile1(Iter first, Iter last)
  220. //-------------------------------------------------------------------
  221. // Returns the first quartile of the elements in [first, last).
  222. //-------------------------------------------------------------------
  223. {
  224.     return quantile(first, last, .25);
  225. };
  226.  
  227. template<class Iter>
  228. inline double quartile3(Iter first, Iter last)
  229. //-------------------------------------------------------------------
  230. // Returns the third quartile of the elements in [first, last).
  231. //-------------------------------------------------------------------
  232. {
  233.     return quantile(first, last, .75);
  234. };
  235.  
  236. template<class Iter>
  237. inline double range(Iter first, Iter last)
  238. //-------------------------------------------------------------------
  239. // Returns the range of the elements in [first, last).
  240. //-------------------------------------------------------------------
  241. {
  242.     return *std::max_element(first, last) - *std::min_element(first, last);
  243. };
  244.  
  245. template<class Iter>
  246. inline double rms(Iter first, Iter last)
  247. //-------------------------------------------------------------------
  248. // Returns the root mean square of the elements in [first, last).
  249. //-------------------------------------------------------------------
  250. {
  251.     return sqrt(sum2(first, last) / length(first, last));
  252. };
  253.  
  254. template<class Iter>
  255. inline double skewness(Iter first, Iter last)
  256. //-------------------------------------------------------------------
  257. // Returns the skewness of the elements in
  258. // [first, last).
  259. //-------------------------------------------------------------------
  260. {
  261.     return centmom(first, last, 3) / pow(stddev(first, last), 3);
  262. };
  263.  
  264. template<class Iter>
  265. inline double stddev(Iter first, Iter last)
  266. //-------------------------------------------------------------------
  267. // Returns the sample standard deviation of the elements in
  268. // [first, last).
  269. //-------------------------------------------------------------------
  270. {
  271.     return sqrt(variance(first, last));
  272. };
  273.  
  274. template<class Iter>
  275. inline double stderrmean(Iter first, Iter last)
  276. //-------------------------------------------------------------------
  277. // Returns the standard error of the mean of the elements in
  278. // [first, last).
  279. //-------------------------------------------------------------------
  280. {
  281.     return sqrt(variance(first, last) / (double)length(first, last));
  282. };
  283.  
  284. template<class Iter>
  285. inline double sum(Iter first, Iter last)
  286. //-------------------------------------------------------------------
  287. // Returns the sum of the elements in [first, last).
  288. //-------------------------------------------------------------------
  289. {
  290.     return std::accumulate(first, last, 0.0);
  291. };
  292.  
  293. template<class Iter>
  294. inline double sum2(Iter first, Iter last)
  295. //-------------------------------------------------------------------
  296. // Returns the sum of the squares of the elements in [first, last).
  297. //-------------------------------------------------------------------
  298. {
  299.     return std::inner_product(first, last, first, (double)0.0);
  300. };
  301.  
  302. template<class Iter1, class Iter2>
  303. inline double sump(Iter1 first1, Iter1 last1, Iter2 first2)
  304. //-------------------------------------------------------------------
  305. // Returns the sum of the products of the corresponding elements in
  306. // [first1, last1) and [first2, first2 + (last1 - first1) ).
  307. //-------------------------------------------------------------------
  308. {
  309.     return std::inner_product(first1, last1, first2, (double)0.0);
  310. };
  311.  
  312. template<class Iter1, class Iter2>
  313. double sumproduct(Iter1 first1, Iter1 last1, Iter2 first2)
  314. //-------------------------------------------------------------------
  315. // Returns the sums of products of the elements in [first, last).
  316. //-------------------------------------------------------------------
  317. {
  318.     int n = length(first1, last1);
  319.     Vector<double> r1(n);
  320.     std::copy(first1, last1, r1.begin());
  321.     Vector<double> r2(n);
  322.     std::copy(first2, first2 + n, r2.begin());
  323.     
  324.     resid(r1.begin(), r1.end(), mean(r1.begin(), r1.end()));
  325.     resid(r2.begin(), r2.end(), mean(r2.begin(), r2.end()));
  326.     
  327.     return std::inner_product(r1.begin(), r1.end(), r2.begin(), 0.0);
  328. };
  329.  
  330. template<class Iter>
  331. double sumsquare(Iter first, Iter last)
  332. //-------------------------------------------------------------------
  333. // Returns the sums of squares of the elements in [first, last).
  334. //-------------------------------------------------------------------
  335. {
  336.     Vector<double> r(length(first, last));
  337.     std::copy(first, last, r.begin());
  338.  
  339.     resid(r.begin(), r.end(), mean(first, last));
  340.     
  341.     double ep = std::accumulate(r.begin(), r.end(), 0.0);
  342.     double s2 = sum2(r.begin(), r.end());
  343.     
  344.     return s2 - ep * ep / (double)r.size();
  345. };
  346.  
  347. template<class Iter>
  348. inline double trimmean(Iter first, Iter last, double f)
  349. //-------------------------------------------------------------------
  350. // Returns the mean of the elements in [first, last) with f*100
  351. // percent of the lowest and highest terms removed.
  352. //-------------------------------------------------------------------
  353. {
  354.     return trimmean(first, last, f, f);
  355. };
  356.  
  357. template<class Iter>
  358. double trimmean(Iter first, Iter last, double f1, double f2)
  359. //-------------------------------------------------------------------
  360. // Returns the mean of the elements in [first, last) with f1*100
  361. // percent of the lowest terms and f2*100 percent of the highest
  362. // terms removed.
  363. //-------------------------------------------------------------------
  364. {
  365.     int n = length(first, last);
  366.     Vector<double> v(n);
  367.     std::copy(first, last, v.begin());
  368.     int n1 = floor(f1 * n), n2 = floor(f2 * n);
  369.     
  370.     std::sort(v.begin(), v.end());
  371.     
  372.     return mean(v.begin() + n1, v.end() - n2);
  373. };
  374.  
  375. template<class Iter>
  376. inline double variance(Iter first, Iter last)
  377. //-------------------------------------------------------------------
  378. // Returns the sample variance of the elements in [first, last).
  379. //-------------------------------------------------------------------
  380. {
  381.     return sumsquare(first, last) / (double)(length(first, last) - 1);
  382. };
  383.  
  384. NUM_END
  385.  
  386. #endif
  387.  
  388. //===================================================================
  389. // Revision History
  390. //
  391. // Version 1.0 - 08/28/1998 - New.
  392. // Version 1.1 - 04/10/1999 - Added std scope resolution.
  393. //             - 04/10/1999 - Added Numerics namespace.
  394. //             - 04/10/1999 - Made length more robust.
  395. //             - 04/10/1999 - Change implementation of mean
  396. //===================================================================
  397.