X-Git-Url: https://git.creatis.insa-lyon.fr/pubgit/?p=CreaPhase.git;a=blobdiff_plain;f=octave_packages%2Fm%2Fstatistics%2Fbase%2Fmedian.m;fp=octave_packages%2Fm%2Fstatistics%2Fbase%2Fmedian.m;h=ab6aa656c7692fb7e36428aca4ace3cf50b4ed0d;hp=0000000000000000000000000000000000000000;hb=1c0469ada9531828709108a4882a751d2816994a;hpb=63de9f36673d49121015e3695f2c336ea92bc278 diff --git a/octave_packages/m/statistics/base/median.m b/octave_packages/m/statistics/base/median.m new file mode 100644 index 0000000..ab6aa65 --- /dev/null +++ b/octave_packages/m/statistics/base/median.m @@ -0,0 +1,123 @@ +## Copyright (C) 1996-2012 John W. Eaton +## Copyright (C) 2009-2010 VZLU Prague +## +## This file is part of Octave. +## +## Octave is free software; you can redistribute it and/or modify it +## under the terms of the GNU General Public License as published by +## the Free Software Foundation; either version 3 of the License, or (at +## your option) any later version. +## +## Octave is distributed in the hope that it will be useful, but +## WITHOUT ANY WARRANTY; without even the implied warranty of +## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU +## General Public License for more details. +## +## You should have received a copy of the GNU General Public License +## along with Octave; see the file COPYING. If not, see +## . + +## -*- texinfo -*- +## @deftypefn {Function File} {} median (@var{x}) +## @deftypefnx {Function File} {} median (@var{x}, @var{dim}) +## Compute the median value of the elements of the vector @var{x}. +## If the elements of @var{x} are sorted, the median is defined +## as +## @tex +## $$ +## {\rm median} (x) = +## \cases{x(\lceil N/2\rceil), & $N$ odd;\cr +## (x(N/2)+x(N/2+1))/2, & $N$ even.} +## $$ +## @end tex +## @ifnottex +## +## @example +## @group +## x(ceil(N/2)) N odd +## median (x) = +## (x(N/2) + x((N/2)+1))/2 N even +## @end group +## @end example +## +## @end ifnottex +## If @var{x} is a matrix, compute the median value for each +## column and return them in a row vector. If the optional @var{dim} +## argument is given, operate along this dimension. +## @seealso{mean, mode} +## @end deftypefn + +## Author: jwe + +function retval = median (x, dim) + + if (nargin != 1 && nargin != 2) + print_usage (); + endif + + if (! (isnumeric (x) || islogical (x))) + error ("median: X must be a numeric vector or matrix"); + endif + + if (isempty (x)) + error ("median: X cannot be an empty matrix"); + endif + + nd = ndims (x); + sz = size (x); + if (nargin < 2) + ## Find the first non-singleton dimension. + (dim = find (sz > 1, 1)) || (dim = 1); + else + if (!(isscalar (dim) && dim == fix (dim)) + || !(1 <= dim && dim <= nd)) + error ("median: DIM must be an integer and a valid dimension"); + endif + endif + + n = sz(dim); + k = floor ((n+1) / 2); + if (mod (n, 2) == 1) + retval = nth_element (x, k, dim); + else + retval = mean (nth_element (x, k:k+1, dim), dim); + endif + ## Inject NaNs where needed, to be consistent with Matlab. + retval(any (isnan (x), dim)) = NaN; + +endfunction + + +%!test +%! x = [1, 2, 3, 4, 5, 6]; +%! x2 = x'; +%! y = [1, 2, 3, 4, 5, 6, 7]; +%! y2 = y'; +%! +%! assert(median (x) == median (x2) && median (x) == 3.5); +%! assert(median (y) == median (y2) && median (y) == 4); +%! assert(median ([x2, 2*x2]) == [3.5, 7]); +%! assert(median ([y2, 3*y2]) == [4, 12]); + +%% Test multidimensional arrays (bug #35679) +%!shared a, b, x, y +%! rand ("seed", 2); +%! a = rand (2,3,4,5); +%! b = rand (3,4,6,5); +%! x = sort (a, 4); +%! y = sort (b, 3); +%!assert (median (a, 4), x(:, :, :, 3)); +%!assert (median (b, 3), (y(:, :, 3, :) + y(:, :, 4, :))/2); + +%!assert(median (single([1,2,3])), single(2)); +%!assert(median ([1,2,NaN;4,5,6;NaN,8,9]), [NaN, 5, NaN]); + +%% Test input validation +%!error median (); +%!error median (1, 2, 3); +%!error median ({1:5}); +%!error median (['A'; 'B']); +%!error median (1, ones(2,2)); +%!error median (1, 1.5); +%!error median (1, 0); +