1## Copyright (C) 2012 Rik Wehbring 2## Copyright (C) 2007-2016 David Bateman 3## 4## This program is free software: you can redistribute it and/or 5## modify it under the terms of the GNU General Public License as 6## published by the Free Software Foundation, either version 3 of the 7## License, or (at your option) any later version. 8## 9## This program is distributed in the hope that it will be useful, but 10## WITHOUT ANY WARRANTY; without even the implied warranty of 11## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU 12## General Public License for more details. 13## 14## You should have received a copy of the GNU General Public License 15## along with this program; see the file COPYING. If not, see 16## <http://www.gnu.org/licenses/>. 17 18## -*- texinfo -*- 19## @deftypefn {} {} unidpdf (@var{x}, @var{n}) 20## For each element of @var{x}, compute the probability density function (PDF) 21## at @var{x} of a discrete uniform distribution which assumes 22## the integer values 1--@var{n} with equal probability. 23## 24## Warning: The underlying implementation uses the double class and will only 25## be accurate for @var{n} < @code{flintmax} (@w{@math{2^{53}}} on 26## IEEE 754 compatible systems). 27## @end deftypefn 28 29function pdf = unidpdf (x, n) 30 31 if (nargin != 2) 32 print_usage (); 33 endif 34 35 if (! isscalar (n)) 36 [retval, x, n] = common_size (x, n); 37 if (retval > 0) 38 error ("unidpdf: X and N must be of common size or scalars"); 39 endif 40 endif 41 42 if (iscomplex (x) || iscomplex (n)) 43 error ("unidpdf: X and N must not be complex"); 44 endif 45 46 if (isa (x, "single") || isa (n, "single")) 47 pdf = zeros (size (x), "single"); 48 else 49 pdf = zeros (size (x)); 50 endif 51 52 k = isnan (x) | ! (n > 0 & n == fix (n)); 53 pdf(k) = NaN; 54 55 k = ! k & (x >= 1) & (x <= n) & (x == fix (x)); 56 if (isscalar (n)) 57 pdf(k) = 1 / n; 58 else 59 pdf(k) = 1 ./ n(k); 60 endif 61 62endfunction 63 64 65%!shared x,y 66%! x = [-1 0 1 2 10 11]; 67%! y = [0 0 0.1 0.1 0.1 0]; 68%!assert (unidpdf (x, 10*ones (1,6)), y) 69%!assert (unidpdf (x, 10), y) 70%!assert (unidpdf (x, 10*[0 NaN 1 1 1 1]), [NaN NaN y(3:6)]) 71%!assert (unidpdf ([x, NaN], 10), [y, NaN]) 72 73## Test class of input preserved 74%!assert (unidpdf (single ([x, NaN]), 10), single ([y, NaN])) 75%!assert (unidpdf ([x, NaN], single (10)), single ([y, NaN])) 76 77## Test input validation 78%!error unidpdf () 79%!error unidpdf (1) 80%!error unidpdf (1,2,3) 81%!error unidpdf (ones (3), ones (2)) 82%!error unidpdf (ones (2), ones (3)) 83%!error unidpdf (i, 2) 84%!error unidpdf (2, i) 85