Mercurial > hg > octave-lyh

comparison scripts/statistics/base/std.m @ 11436:e151e23f73bc

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
Overhaul base statistics functions and documentation of same. Add or improve input validation. Add input validation tests. Add functional tests. Improve or re-write documentation strings.
author Rik <octave@nomad.inbox5.com>
date Mon, 03 Jan 2011 21:23:08 -0800
parents 693e22af08ae
children fd0a3ac60b0e
comparison
equal deleted inserted replaced
11435:20f53b3a558f 11436:e151e23f73bc
19 19
20 ## -*- texinfo -*- 20 ## -*- texinfo -*-
21 ## @deftypefn {Function File} {} std (@var{x}) 21 ## @deftypefn {Function File} {} std (@var{x})
22 ## @deftypefnx {Function File} {} std (@var{x}, @var{opt}) 22 ## @deftypefnx {Function File} {} std (@var{x}, @var{opt})
23 ## @deftypefnx {Function File} {} std (@var{x}, @var{opt}, @var{dim}) 23 ## @deftypefnx {Function File} {} std (@var{x}, @var{opt}, @var{dim})
24 ## If @var{x} is a vector, compute the standard deviation of the elements 24 ## Compute the standard deviation of the elements of the vector @var{x}.
25 ## of @var{x}.
26 ## @tex 25 ## @tex
27 ## $$ 26 ## $$
28 ## {\rm std} (x) = \sigma (x) = \sqrt{{\sum_{i=1}^N (x_i - \bar{x})^2 \over N - 1}} 27 ## {\rm std} (x) = \sigma = \sqrt{{\sum_{i=1}^N (x_i - \bar{x})^2 \over N - 1}}
29 ## $$ 28 ## $$
30 ## where $\bar{x}$ is the mean value of $x$. 29 ## where $\bar{x}$ is the mean value of $x$ and $N$ is the number of elements.
31 ## @end tex 30 ## @end tex
32 ## @ifnottex 31 ## @ifnottex
33 ## 32 ##
34 ## @example 33 ## @example
35 ## @group 34 ## @group
36 ## std (x) = sqrt (sumsq (x - mean (x)) / (n - 1)) 35 ## std (x) = sqrt ( 1/(N-1) SUM_i (x(i) - mean(x))^2 )
37 ## @end group 36 ## @end group
38 ## @end example 37 ## @end example
39 ## 38 ##
39 ## @noindent
40 ## where @math{N} is the number of elements.
40 ## @end ifnottex 41 ## @end ifnottex
42 ##
41 ## If @var{x} is a matrix, compute the standard deviation for 43 ## If @var{x} is a matrix, compute the standard deviation for
42 ## each column and return them in a row vector. 44 ## each column and return them in a row vector.
43 ## 45 ##
44 ## The argument @var{opt} determines the type of normalization to use. Valid 46 ## The argument @var{opt} determines the type of normalization to use.
45 ## values are 47 ## Valid values are
46 ## 48 ##
47 ## @table @asis 49 ## @table @asis
48 ## @item 0: 50 ## @item 0:
49 ## normalizes with @math{N-1}, provides the square root of best unbiased 51 ## normalize with @math{N-1}, provides the square root of the best unbiased
50 ## estimator of the variance [default] 52 ## estimator of the variance [default]
51 ## 53 ##
52 ## @item 1: 54 ## @item 1:
53 ## normalizes with @math{N}, this provides the square root of the second 55 ## normalize with @math{N}, this provides the square root of the second
54 ## moment around the mean 56 ## moment around the mean
55 ## @end table 57 ## @end table
56 ## 58 ##
57 ## The third argument @var{dim} determines the dimension along which the 59 ## If the optional argument @var{dim} is given, operate along this dimension.
58 ## standard 60 ## @seealso{var, range, iqr, mean, median}
59 ## deviation is calculated.
60 ## @seealso{mean, median}
61 ## @end deftypefn 61 ## @end deftypefn
62 62
63 ## Author: jwe 63 ## Author: jwe
64 64
65 function retval = std (a, opt, dim) 65 function retval = std (x, opt = 0, dim)
66 66
67 if (nargin < 1 || nargin > 3) 67 if (nargin < 1 || nargin > 3)
68 print_usage (); 68 print_usage ();
69 endif 69 endif
70
71 if (! (isnumeric (x)))
72 error ("std: X must be a numeric vector or matrix");
73 endif
74
75 if (isempty (opt))
76 opt = 0;
77 endif
78 if (opt != 0 && opt != 1)
79 error ("std: normalization OPT must be 0 or 1");
80 endif
81
82 sz = size (x);
70 if (nargin < 3) 83 if (nargin < 3)
71 dim = find (size (a) > 1, 1); 84 ## Find the first non-singleton dimension.
85 dim = find (sz > 1, 1);
72 if (isempty (dim)) 86 if (isempty (dim))
73 dim = 1; 87 dim = 1;
74 endif 88 endif
75 endif 89 endif
76 if (nargin < 2 || isempty (opt))
77 opt = 0;
78 endif
79 90
80 n = size (a, dim); 91 n = size (x, dim);
81 if (n == 1) 92 if (n == 1)
82 retval = zeros (size (a)); 93 retval = zeros (sz);
83 elseif (numel (a) > 0) 94 elseif (numel (x) > 0)
84 retval = sqrt (sumsq (center (a, dim), dim) / (n + opt - 1)); 95 retval = sqrt (sumsq (center (x, dim), dim) / (n - 1 + opt));
85 else 96 else
86 error ("std: x must not be empty"); 97 error ("std: X must not be empty");
87 endif 98 endif
88 99
89 endfunction 100 endfunction
101
90 102
91 %!test 103 %!test
92 %! x = ones (10, 2); 104 %! x = ones (10, 2);
93 %! y = [1, 3]; 105 %! y = [1, 3];
94 %! assert(std (x) == [0, 0] && abs (std (y) - sqrt (2)) < sqrt (eps)); 106 %! assert(std (x) == [0, 0] && abs (std (y) - sqrt (2)) < sqrt (eps));
95 %! assert (std (x, 0, 3), zeros (10, 2)) 107 %! assert (std (x, 0, 3), zeros (10, 2))
96 %! assert (std (ones (3, 1, 2), 0, 2), zeros (3, 1, 2)) 108 %! assert (std (ones (3, 1, 2), 0, 2), zeros (3, 1, 2))
97 109
110 %% Test input validation
98 %!error std (); 111 %!error std ();
112 %!error std (1, 2, 3, 4);
113 %!error std (true(1,2))
114 %!error std (1, -1);
115 %!error std ([], 1);
99 116
100 %!error std (1, 2, 3, 4);