Mercurial > octave-nkf

comparison scripts/statistics/base/moment.m @ 17873:58b39152b0f6

Find changesets by keywords (author, files, the commit message), revision number or hash, or revset expression.
moment.m: Compute central, rather than raw, moment for ML compatibility (bug #36718). * NEWS: Announce change in definition of moment. * scripts/statistics/base/moment.m: Default to calculating central moment. Add "r" option for calculating raw moment.
author Rik <rik@octave.org>
date Thu, 07 Nov 2013 09:39:49 -0800
parents d63878346099
children ed2ef5d96929
comparison
equal deleted inserted replaced
17872:7d9a4eef8022 17873:58b39152b0f6
20 ## @deftypefn {Function File} {} moment (@var{x}, @var{p}) 20 ## @deftypefn {Function File} {} moment (@var{x}, @var{p})
21 ## @deftypefnx {Function File} {} moment (@var{x}, @var{p}, @var{type}) 21 ## @deftypefnx {Function File} {} moment (@var{x}, @var{p}, @var{type})
22 ## @deftypefnx {Function File} {} moment (@var{x}, @var{p}, @var{dim}) 22 ## @deftypefnx {Function File} {} moment (@var{x}, @var{p}, @var{dim})
23 ## @deftypefnx {Function File} {} moment (@var{x}, @var{p}, @var{type}, @var{dim}) 23 ## @deftypefnx {Function File} {} moment (@var{x}, @var{p}, @var{type}, @var{dim})
24 ## @deftypefnx {Function File} {} moment (@var{x}, @var{p}, @var{dim}, @var{type}) 24 ## @deftypefnx {Function File} {} moment (@var{x}, @var{p}, @var{dim}, @var{type})
25 ## Compute the @var{p}-th moment of the vector @var{x} about zero. 25 ## Compute the @var{p}-th central moment of the vector @var{x}.
26 ## @tex 26 ##
27 ## $$ 27 ## @tex
28 ## {\rm moment} (x) = { \sum_{i=1}^N {x_i}^p \over N } 28 ## $$
29 ## $$ 29 ## {\sum_{i=1}^N (x_i - \bar{x})^p \over N}
30 ## @end tex 30 ## $$
31 ## @ifnottex 31 ## @end tex
32 ## 32 ## @ifnottex
33 ## @example 33 ##
34 ## @group 34 ## @example
35 ## moment (x) = 1/N SUM_i x(i)^p 35 ## @group
36 ## @end group 36 ## 1/N SUM_i (x(i) - mean(x))^p
37 ## @end example 37 ## @end group
38 ## 38 ## @end example
39 ## @end ifnottex 39 ##
40 ## 40 ## @end ifnottex
41 ## If @var{x} is a matrix, return the row vector containing the 41 ##
42 ## @var{p}-th moment of each column. 42 ## If @var{x} is a matrix, return the row vector containing the @var{p}-th
43 ## central moment of each column.
43 ## 44 ##
44 ## The optional string @var{type} specifies the type of moment to be computed. 45 ## The optional string @var{type} specifies the type of moment to be computed.
45 ## Valid options are: 46 ## Valid options are:
46 ## 47 ##
47 ## @table @asis 48 ## @table @asis
48 ## @item @qcode{"c"} 49 ## @item @qcode{"c"}
49 ## Central Moment. The moment about the mean defined as 50 ## Central Moment (default).
50 ## @tex 51 ##
51 ## $$ 52 ## @item @qcode{"a"}
52 ## {\sum_{i=1}^N (x_i - \bar{x})^p \over N} 53 ## @itemx @qcode{"ac"}
53 ## $$ 54 ## Absolute Central Moment. The moment about the mean ignoring sign
54 ## @end tex 55 ## defined as
55 ## @ifnottex 56 ## @tex
56 ## 57 ## $$
57 ## @example 58 ## {\sum_{i=1}^N {\left| x_i - \bar{x} \right|}^p \over N}
58 ## @group 59 ## $$
59 ## 1/N SUM_i (x(i) - mean(x))^p 60 ## @end tex
60 ## @end group 61 ## @ifnottex
61 ## @end example 62 ##
62 ## 63 ## @example
63 ## @end ifnottex 64 ## @group
64 ## 65 ## 1/N SUM_i (abs (x(i) - mean(x)))^p
65 ## @item @qcode{"a"} 66 ## @end group
66 ## Absolute Moment. The moment about zero ignoring sign defined as 67 ## @end example
68 ##
69 ## @end ifnottex
70 ##
71 ## @item @qcode{"r"}
72 ## Raw Moment. The moment about zero defined as
73 ##
74 ## @tex
75 ## $$
76 ## {\rm moment} (x) = { \sum_{i=1}^N {x_i}^p \over N }
77 ## $$
78 ## @end tex
79 ## @ifnottex
80 ##
81 ## @example
82 ## @group
83 ## moment (x) = 1/N SUM_i x(i)^p
84 ## @end group
85 ## @end example
86 ##
87 ## @end ifnottex
88 ##
89 ## @item @qcode{"ar"}
90 ## Absolute Raw Moment. The moment about zero ignoring sign defined as
67 ## @tex 91 ## @tex
68 ## $$ 92 ## $$
69 ## {\sum_{i=1}^N {\left| x_i \right|}^p \over N} 93 ## {\sum_{i=1}^N {\left| x_i \right|}^p \over N}
70 ## $$ 94 ## $$
71 ## @end tex 95 ## @end tex
76 ## 1/N SUM_i ( abs (x(i)) )^p 100 ## 1/N SUM_i ( abs (x(i)) )^p
77 ## @end group 101 ## @end group
78 ## @end example 102 ## @end example
79 ## 103 ##
80 ## @end ifnottex 104 ## @end ifnottex
81 ##
82 ## @item @qcode{"ac"}
83 ## Absolute Central Moment. Defined as
84 ## @tex
85 ## $$
86 ## {\sum_{i=1}^N {\left| x_i - \bar{x} \right|}^p \over N}
87 ## $$
88 ## @end tex
89 ## @ifnottex
90 ##
91 ## @example
92 ## @group
93 ## 1/N SUM_i ( abs (x(i) - mean(x)) )^p
94 ## @end group
95 ## @end example
96 ##
97 ## @end ifnottex
98 ## @end table 105 ## @end table
99 ##
100 ## If the optional argument @var{dim} is given, operate along this dimension. 106 ## If the optional argument @var{dim} is given, operate along this dimension.
101 ## 107 ##
102 ## If both @var{type} and @var{dim} are given they may appear in any order. 108 ## If both @var{type} and @var{dim} are given they may appear in any order.
103 ## @seealso{var, skewness, kurtosis} 109 ## @seealso{var, skewness, kurtosis}
104 ## @end deftypefn 110 ## @end deftypefn
113 119
114 if (nargin < 2 || nargin > 4) 120 if (nargin < 2 || nargin > 4)
115 print_usage (); 121 print_usage ();
116 endif 122 endif
117 123
118 if (!(isnumeric (x) || islogical (x)) || isempty (x)) 124 if (! (isnumeric (x) || islogical (x)) || isempty (x))
119 error ("moment: X must be a non-empty numeric matrix or vector"); 125 error ("moment: X must be a non-empty numeric matrix or vector");
120 endif 126 endif
121 127
122 if (! (isnumeric (p) && isscalar (p))) 128 if (! (isnumeric (p) && isscalar (p)))
123 error ("moment: P must be a numeric scalar"); 129 error ("moment: P must be a numeric scalar");
160 endif 166 endif
161 endif 167 endif
162 168
163 n = sz(dim); 169 n = sz(dim);
164 170
165 if (any (type == "c")) 171 if (! any (type == "r"))
166 x = center (x, dim); 172 x = center (x, dim);
167 endif 173 endif
168 if (any (type == "a")) 174 if (any (type == "a"))
169 x = abs (x); 175 x = abs (x);
170 endif 176 endif
174 endfunction 180 endfunction
175 181
176 182
177 %!test 183 %!test
178 %! x = rand (10); 184 %! x = rand (10);
179 %! assert (moment (x,1), mean (x), 1e1*eps); 185 %! assert (moment (x,1), mean (center (x)));
180 %! assert (moment (x,2), meansq (x), 1e1*eps); 186 %! assert (moment (x,2), meansq (center (x)));
181 %! assert (moment (x,1,2), mean (x,2), 1e1*eps); 187 %! assert (moment (x,1,2), mean (center (x, 2), 2));
182 %! assert (moment (x,1,"c"), mean (center (x)), 1e1*eps); 188 %! assert (moment (x,1,"a"), mean (abs (center (x))));
183 %! assert (moment (x,1,"a"), mean (abs (x)), 1e1*eps); 189 %! assert (moment (x,1,"r"), mean (x));
184 190 %! assert (moment (x,1,"ar"), mean (abs (x)));
185 %!assert (moment (single ([1 2 3]), 1), single (2)) 191
192 %!assert (moment (single ([1 2 3]), 1, "r"), single (2))
186 193
187 %% Test input validation 194 %% Test input validation
188 %!error moment () 195 %!error moment ()
189 %!error moment (1) 196 %!error moment (1)
190 %!error moment (1, 2, 3, 4, 5) 197 %!error moment (1, 2, 3, 4, 5)
191 %!error moment (['A'; 'B'], 2) 198 %!error <X must be a non-empty numeric matrix> moment (['A'; 'B'], 2)
192 %!error moment (ones (2,0,3), 2) 199 %!error <X must be a non-empty numeric matrix> moment (ones (2,0,3), 2)
193 %!error moment (1, true) 200 %!error <P must be a numeric scalar> moment (1, true)
194 %!error moment (1, ones (2,2)) 201 %!error <P must be a numeric scalar> moment (1, ones (2,2))
195 %!error moment (1, 2, 3, 4) 202 %!error <expecting TYPE to be a string> moment (1, 2, 3, 4)
196 %!error moment (1, 2, ones (2,2)) 203 %!error <DIM must be an integer and a valid dimension> moment (1, 2, ones (2,2))
197 %!error moment (1, 2, 1.5) 204 %!error <DIM must be an integer and a valid dimension> moment (1, 2, 1.5)
198 %!error moment (1, 2, 4) 205 %!error <DIM must be an integer and a valid dimension> moment (1, 2, 4)
199 206