annotate scripts/specfun/betaincinv.m @ 24927:c280560d9c96 stable

Overhaul special functions modified by GSOC2018 project. * NEWS: Add note about new functions added. Add note explaining the changes done to existing functions. * __betainc__.cc: Renamed from __betainc_lentz__.cc. Use standard GPL v3 copyright block. Add missing #include "dNDArray.h". Add one-line texinfo documentation for internal function. Remove fourth input argument to function. Use is_single_type() to decide whether to operate with FloatNDArray or NDArray. Delete temporary variables x_arg_s, a_arg_s, b_arg_s used in input validation. Rename len_x, len_a, len_b to numel_[xab] for clarity. Remove input validation that numel match (internal function, no need). Rename function output to retval according to Octave coding conventions. Use more spacing and newlines for readability of code. * __expint__.cc: Renamed from __expint_lentz__.cc. Use standard GPL v3 copyright block. Cull list of #includes to just the necessary ones. Use Complex and FloatComplex typedefs defined by Octave. Add one-line texinfo documentation for internal function. Remove second input argument to function. Use is_single_type() to decide whether to operate with FloatComplexNDArray or ComplexNDArray. Delete temporary variables x_arg_s used in input validation. Rename len_x to numel_x for clarity. Use constructor with dim_vector and scalar value rather than fill() after creating array. Rename function output to retval according to Octave coding conventions. Use more spacing and newlines for readability of code. * __gammainc__.cc: Renamed from __gammainc_lentz__.cc. Use standard GPL v3 copyright block. Add one-line texinfo documentation for internal function. Remove third input argument to function. Remove input validation that numel match (internal function, no need). Use is_single_type() to decide whether to operate with FloatNDArray or NDArray. Delete temporary variables x_arg, a_arg used in input validation. Use constructor with dim_vector and scalar value rather than fill() after creating array. Rename function output to retval according to Octave coding conventions. Use more spacing and newlines for readability of code. * __betainc_lentz__.cc, __expint_lentz__.cc, __gammainc_lentz__.cc: Removed. * libinterp/corefcn/module.mk: Add renamed functions to build system. * betainc.m: Use Octave standard GPL block. Rewrite parts of docstring. Don't use array brackets around single output of function. Remove isscalar checks on inputs because common_size() function will already handle it. Use capital variable names in error messages to match documentation as displayed in terminal. Reshape all inputs in to column vectors quickly so that input validation tests that depend on all/any will pass with N-D arrays. Add comments to code. Check for specific error messages in input validation BIST tests. * betaincinv.m: Use Octave standard GPL block. Rewrite parts of docstring. Don't use array brackets around single output of function. Remove isscalar checks on inputs because common_size() function will already handle it. Use capital variable names in error messages to match documentation as displayed in terminal. Reshape all inputs in to column vectors quickly so that input validation tests that depend on all/any will pass with N-D arrays. Put most common case of tail ("lower") first in if/elseif trees. Call functions directly with function handle rather than using unnecessary feval() call. Use numel in preference to length. Rename variable i_miss to todo for clarity. Add comments to code. Check for specific error messages in input validation BIST tests. * cosint.m: Use Octave standard GPL block. Rewrite parts of docstring. Don't use array brackets around single output of function. Remove isscalar checks on inputs because common_size() function will already handle it. Add input validation check for isnumeric value. Convert integer classes to double before proceeding. Rename i_miss to todo for clarity. Use isinf to detect both -Inf and +Inf rather than separate tests. Use ++it in while loop conditional to shorten loop blocks. Add Input validation BIST tests. * expint.m: Remove Sylvain who was not actually an author on this file. Rewrite parts of docstring. Rename variable sparse_x to orig_sparse. Eliminate temporary variables res_tmp, x_s_tmp, ssum_tmp. Rename i_miss to todo. Use Octave coding conventions throughout. Add comments to code. * gammainc.m: Use Octave standard GPL block. Rewrite parts of docstring. Remove isscalar checks on inputs because common_size() function will already handle it. Use capital variable names in error messages to match documentation as displayed in terminal. Reshape all inputs in to column vectors quickly so that input validation tests that depend on all/any will pass with N-D arrays. Put most common case of tail ("lower") first in if/elseif trees. Add input validation of tail. Rename variable ii to idx for clarity. Rename variable i_done to todo and switch polarity so that the '!' operator is not required every time the variable is updated. Use indexing and direct assignment to update todo rather than logical operator '&' which is slower. Use tolower on tail variable and then switch strcmpi calls to strcmp. Reformat %!test blocks in to %!assert blocks to be more compact. Check for specific error messages in input validation BIST tests. * gammaincinv.m: Use Octave standard GPL block. Rewrite parts of docstring. Don't use array brackets around single output of function. Remove isscalar checks on inputs because common_size() function will already handle it. Add input validation check for iscomplex value. Use capital variable names in error messages to match documentation as displayed in terminal. Reshape all inputs in to column vectors quickly so that input validation tests that depend on all/any will pass with N-D arrays. Rename i_miss to todo. Use numel in preference to length. Call functions directly with function handle rather than using unnecessary feval() call.Rename variable i_miss to todo for clarity. Use it++ in while loop conditional to shorten loop blocks. Add comments to code. Check for specific error messages in input validation BIST tests. Add input validation BIST tests for all error messages. * sinint.m: Use Octave standard GPL block. Rewrite parts of docstring. Add input validation for isnumeric. Convert integers to double for calculation. Reshape input to column vector. Rename variable sz to orig_sz for clarity. rename i_miss to todo. Reformat BIST tests to mak them more compact. Add input validation BIST tests.
author Rik <rik@octave.org>
date Mon, 19 Mar 2018 10:01:48 -0700
parents bd89440407aa
children f80e68529bbf
Ignore whitespace changes - Everywhere: Within whitespace: At end of lines:
rev   line source
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
1 ## Copyright (C) 2017 Michele Ginesi
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
2 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
3 ## This file is part of Octave.
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
4 ##
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
5 ## Octave is free software: you can redistribute it and/or modify it
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
6 ## under the terms of the GNU General Public License as published by
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
7 ## the Free Software Foundation, either version 3 of the License, or
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
8 ## (at your option) any later version.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
9 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
10 ## Octave is distributed in the hope that it will be useful, but
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
11 ## WITHOUT ANY WARRANTY; without even the implied warranty of
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
13 ## GNU General Public License for more details.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
14 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
15 ## You should have received a copy of the GNU General Public License
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
16 ## along with Octave; see the file COPYING. If not, see
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
17 ## <https://www.gnu.org/licenses/>.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
18
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
19 ## Author: Michele Ginesi <michele.ginesi@gmail.com>
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
20
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
21 ## -*- texinfo -*-
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
22 ## @deftypefn {} {} betaincinv (@var{y}, @var{a}, @var{b})
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
23 ## @deftypefnx {} {} betaincinv (@var{y}, @var{a}, @var{b}, "lower")
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
24 ## @deftypefnx {} {} betaincinv (@var{y}, @var{a}, @var{b}, "upper")
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
25 ## Compute the inverse of the normalized incomplete beta function.
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
26 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
27 ## The normalized incomplete beta function is defined as
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
28 ## @tex
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
29 ## $$
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
30 ## I_x (a, b) = {1 \over {B(a,b)}} \displaystyle{\int_0^x t^{a-1} (1-t)^{b-1} dt}
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
31 ## $$
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
32 ## @end tex
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
33 ## @ifnottex
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
34 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
35 ## @example
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
36 ## @group
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
37 ## x
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
38 ## /
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
39 ## |
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
40 ## I_x (a, b) = | t^(a-1) (1-t)^(b-1) dt
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
41 ## |
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
42 ## /
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
43 ## 0
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
44 ## @end group
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
45 ## @end example
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
46 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
47 ## @end ifnottex
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
48 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
49 ## If two inputs are scalar, then @code{betaincinv (@var{y}, @var{a}, @var{b})}
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
50 ## is returned for each of the other inputs.
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
51 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
52 ## If two or more inputs are not scalar, the sizes of them must agree, and
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
53 ## @code{betaincinv} is applied element-by-element.
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
54 ##
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
55 ## The variable @var{y} must be in the interval [0,1], while @var{a} and
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
56 ## @var{b} must be real and strictly positive.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
57 ##
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
58 ## By default, @var{tail} is @qcode{"lower")} and the inverse of the incomplete
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
59 ## beta function integrated from 0 to @var{x} is computed. If @var{tail} is
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
60 ## @qcode{"upper"} then the complementary function integrated from @var{x} to 1
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
61 ## is inverted.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
62 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
63 ## The function is computed by standard Newton's method, by solving
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
64 ## @tex
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
65 ## $$
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
66 ## y - I_x (a, b) = 0
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
67 ## $$
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
68 ## @end tex
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
69 ## @ifnottex
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
70 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
71 ## @example
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
72 ## @var{y} - betainc (@var{x}, @var{a}, @var{b}) = 0
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
73 ## @end example
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
74 ##
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
75 ## @end ifnottex
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
76 ##
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
77 ## @seealso{betainc, beta, betaln}
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
78 ## @end deftypefn
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
79
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
80 function x = betaincinv (y, a, b, tail = "lower")
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
81
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
82 if (nargin < 3 || nargin > 4)
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
83 print_usage ();
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
84 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
85
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
86 [err, y, a, b] = common_size (y, a, b);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
87 if (err > 0)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
88 error ("betaincinv: Y, A, and B must be of common size or scalars");
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
89 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
90
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
91 if (iscomplex (y) || iscomplex (a) || iscomplex (b))
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
92 error ("betaincinv: all inputs must be real");
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
93 endif
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
94
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
95 ## FIXME: Should there be isnumeric checking? Right now it accepts char
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
96 ## arrays, but then produces a weird error later on.
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
97
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
98 ## Remember original shape of data, but convert to column vector for calcs.
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
99 orig_sz = size (y);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
100 y = y(:);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
101 a = a(:);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
102 b = b(:);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
103
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
104 if (any ((y < 0) | (y > 1)))
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
105 error ("betaincinv: Y must be in the range [0, 1]");
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
106 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
107
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
108 if (any (a <= 0))
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
109 error ("betaincinv: A must be strictly positive");
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
110 endif
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
111
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
112 if (any (b <= 0))
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
113 error ("betaincinv: B must be strictly positive");
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
114 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
115
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
116 ## If any of the arguments is single then the output should be as well.
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
117 if (strcmp (class (y), "single") || strcmp (class (a), "single")
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
118 || strcmp (class (b), "single"))
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
119 y = single (y);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
120 a = single (a);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
121 b = single (b);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
122 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
123
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
124 ## Convert to floating point if necessary
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
125 if (isinteger (y))
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
126 y = double (y);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
127 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
128 if (isinteger (a))
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
129 a = double (a);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
130 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
131 if (isinteger (b))
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
132 b = double (b);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
133 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
134
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
135 ## Initialize output array
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
136 x = zeros (size (y), class (y));
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
137
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
138 ## Parameters for the Newton method
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
139 maxit = 20;
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
140 tol = eps (class (y));
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
141
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
142 if (strcmpi (tail, "lower"))
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
143 p = y;
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
144 q = 1 - y;
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
145 x(y == 0) = 0;
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
146 x(y == 1) = 1;
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
147 elseif (strcmpi (tail, "upper"))
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
148 p = 1 - y;
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
149 q = y;
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
150 x(y == 0) = 1;
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
151 x(y == 1) = 0;
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
152 else
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
153 error ("betaincinv: invalid value for TAIL")
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
154 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
155
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
156 ## Special values have been already computed.
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
157 todo = (y != 0) & (y != 1);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
158
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
159 ## We will invert the lower version for p < 0.5 and the upper otherwise.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
160 i_low = (p < 0.5);
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
161 i_upp = (! i_low);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
162
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
163 idx = todo & i_low;
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
164 if (any (idx));
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
165 n = nnz (idx);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
166 ## Function and derivative of the lower version.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
167 F = @(x, a, b, y) y - betainc (x, a, b);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
168 JF = @(x, a, b) - real (exp ((a-1) .* log (x) + (b-1) .* log1p (-x) + ...
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
169 gammaln (a+b) - gammaln (a) - gammaln (b)));
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
170
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
171 ## Compute the initial guess with a bisection method of 10 steps.
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
172 x0 = bisection_method (F, zeros (n,1), ones (n,1), ...
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
173 a(i_low), b(i_low), p(i_low), 10);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
174
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
175 ## Use Newton's method to iteratively find solution.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
176 x(i_low) = newton_method (F, JF, x0, a(i_low), b(i_low), p(i_low), ...
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
177 tol, maxit);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
178 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
179
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
180 idx = todo & i_upp;
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
181 if (any (idx));
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
182 n = nnz (idx);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
183 ## Function and derivative of the upper version.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
184 F = @(x, a, b, y) y - betainc (x, a, b, "upper");
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
185 JF = @(x, a, b) real (exp ((a-1) .* log (x) + (b-1) .* log1p (-x) + ...
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
186 gammaln (a+b) - gammaln (a) - gammaln (b)));
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
187
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
188 ## Compute the initial guess with a bisection method of 10 steps.
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
189 x0 = bisection_method (F, zeros (n,1), ones (n,1), ...
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
190 a(i_upp), b(i_upp), q(i_upp), 10);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
191
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
192 ## Use Newton's method to iteratively find solution.
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
193 x(i_upp) = newton_method (F, JF, x0, a(i_upp), b(i_upp), q(i_upp), ...
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
194 tol, maxit);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
195 endif
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
196
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
197 ## Restore original shape
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
198 x = reshape (x, orig_sz);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
199
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
200 endfunction
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
201
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
202
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
203 ## Subfunctions: Bisection and Newton Methods
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
204 function xc = bisection_method (F, xl, xr, a, b, y, maxit)
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
205 F_l = F (xl, a, b, y);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
206 F_r = F (xr, a, b, y);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
207 for it = 1:maxit
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
208 xc = (xl + xr) / 2;
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
209 F_c = F (xc, a, b, y);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
210 flag_l = ((F_c .* F_r) < 0);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
211 flag_r = ((F_c .* F_l) < 0);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
212 flag_c = (F_c == 0);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
213 xl(flag_l) = xc(flag_l);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
214 xr(flag_r) = xc(flag_r);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
215 xl(flag_c) = xr(flag_c) = xc(flag_c);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
216 F_l(flag_l) = F_c(flag_l);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
217 F_r(flag_r) = F_c(flag_r);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
218 F_l(flag_c) = F_r(flag_c) = 0;
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
219 endfor
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
220 endfunction
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
221
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
222 function x = newton_method (F, JF, x0, a, b, y, tol, maxit);
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
223 l = numel (y);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
224 res = -F (x0, a, b, y) ./ JF (x0, a, b);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
225 todo = (abs(res) >= tol * abs (x0));
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
226 x = x0;
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
227 it = 0;
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
228 while (any (todo) && (it < maxit))
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
229 it++;
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
230 x(todo) += res(todo);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
231 res(todo) = -F(x(todo), a(todo), b(todo), y(todo)) ...
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
232 ./ JF (x(todo), a(todo), b(todo));
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
233 todo = (abs(res) >= tol * abs (x));
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
234 endwhile
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
235 x += res;
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
236 endfunction
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
237
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
238
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
239 %!test
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
240 %! x = linspace (0.1, 0.9, 11);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
241 %! a = [2, 3, 4];
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
242 %! [x,a,b] = ndgrid (x,a,a);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
243 %! xx = betaincinv (betainc (x, a, b), a, b);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
244 %! assert (xx, x, 3e-15);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
245
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
246 %!test
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
247 %! x = linspace (0.1, 0.9, 11);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
248 %! a = [2, 3, 4];
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
249 %! [x,a,b] = ndgrid (x,a,a);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
250 %! xx = betaincinv (betainc (x, a, b, "upper"), a, b, "upper");
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
251 %! assert (xx, x, 3e-15);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
252
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
253 %!test
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
254 %! x = linspace (0.1, 0.9, 11);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
255 %! a = [0.1:0.1:1];
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
256 %! [x,a,b] = ndgrid (x,a,a);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
257 %! xx = betaincinv (betainc (x, a, b), a, b);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
258 %! assert (xx, x, 3e-15);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
259
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
260 %!test
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
261 %! x = linspace (0.1, 0.9, 11);
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
262 %! a = [0.1:0.1:1];
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
263 %! [x,a,b] = ndgrid (x,a,a);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
264 %! xx = betaincinv (betainc (x, a, b, "upper"), a, b, "upper");
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
265 %! assert (xx, x, 3e-15);
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
266
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
267 ## Test the conservation of the input class
24907
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
268 %!assert (class (betaincinv (0.5, 1, 1)), "double")
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
269 %!assert (class (betaincinv (single (0.5), 1, 1)), "single")
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
270 %!assert (class (betaincinv (0.5, single (1), 1)), "single")
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
271 %!assert (class (betaincinv (int8 (0), 1, 1)), "double")
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
272 %!assert (class (betaincinv (0.5, int8 (1), 1)), "double")
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
273 %!assert (class (betaincinv (int8 (0), single (1), 1)), "single")
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
274 %!assert (class (betaincinv (single (0.5), int8 (1), 1)), "single")
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
275
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
276 ## Test input validation
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
277 %!error betaincinv ()
bd89440407aa Incomplete beta function moved to a .m file, fixing accuracy and
Michele Ginesi <michele.ginesi@gmail.com>
parents:
diff changeset
278 %!error betaincinv (1)
24927
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
279 %!error betaincinv (1,2)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
280 %!error betaincinv (1,2,3,4,5)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
281 %!error <must be of common size or scalars>
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
282 %! betaincinv (ones (2,2), ones (1,2), 1);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
283 %!error <all inputs must be real> betaincinv (0.5i, 1, 2)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
284 %!error <all inputs must be real> betaincinv (0, 1i, 1)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
285 %!error <all inputs must be real> betaincinv (0, 1, 1i)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
286 %!error <Y must be in the range \[0, 1\]> betaincinv (-0.1,1,1)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
287 %!error <Y must be in the range \[0, 1\]> betaincinv (1.1,1,1)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
288 %!error <Y must be in the range \[0, 1\]>
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
289 %! y = ones (1, 1, 2);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
290 %! y(1,1,2) = -1;
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
291 %! betaincinv (y,1,1);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
292 %!error <A must be strictly positive> betaincinv (0.5,0,1)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
293 %!error <A must be strictly positive>
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
294 %! a = ones (1, 1, 2);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
295 %! a(1,1,2) = 0;
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
296 %! betaincinv (1,a,1);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
297 %!error <B must be strictly positive> betaincinv (0.5,1,0)
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
298 %!error <B must be strictly positive>
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
299 %! b = ones (1, 1, 2);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
300 %! b(1,1,2) = 0;
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
301 %! betaincinv (1,1,b);
c280560d9c96 Overhaul special functions modified by GSOC2018 project.
Rik <rik@octave.org>
parents: 24907
diff changeset
302 %!error <invalid value for TAIL> betaincinv (1,2,3, "foobar")