本文整理汇总了Java中org.renjin.sexp.DoubleVector.isFinite方法的典型用法代码示例。如果您正苦于以下问题:Java DoubleVector.isFinite方法的具体用法?Java DoubleVector.isFinite怎么用?Java DoubleVector.isFinite使用的例子?那么恭喜您, 这里精选的方法代码示例或许可以为您提供帮助。您也可以进一步了解该方法所在类org.renjin.sexp.DoubleVector
的用法示例。
在下文中一共展示了DoubleVector.isFinite方法的15个代码示例,这些例子默认根据受欢迎程度排序。您可以为喜欢或者感觉有用的代码点赞,您的评价将有助于系统推荐出更棒的Java代码示例。
示例1: dpois_raw
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double dpois_raw(double x, double lambda, boolean give_log) {
/* x >= 0 ; integer for dpois(), but not e.g. for pgamma()!
lambda >= 0
*/
if (lambda == 0) {
return ((x == 0) ? SignRank.R_D__1(true, give_log) : SignRank.R_D__0(true, give_log));
}
if (!DoubleVector.isFinite(lambda)) {
return SignRank.R_D__0(true, give_log);
}
if (x < 0) {
return (SignRank.R_D__0(true, give_log));
}
if (x <= lambda * Double.MIN_VALUE) {
return (SignRank.R_D_exp(-lambda, true, give_log));
}
if (lambda < x * Double.MIN_VALUE) {
return (SignRank.R_D_exp(-lambda + x * Math.log(lambda) - Gamma.logGamma(x + 1), true, give_log));
}
return (SignRank.R_D_fexp(2 * Math.PI * x, -Binom.stirlerr(x) - Binom.bd0(x, lambda), true, give_log));
}
示例2: bd0
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double bd0(double x, double np) {
double ej, s, s1, v;
int j;
if (!DoubleVector.isFinite(x) || !DoubleVector.isFinite(np) || np == 0.0) {
return DoubleVector.NaN;
}
if (Math.abs(x - np) < 0.1 * (x + np)) {
v = (x - np) / (x + np);
s = (x - np) * v;/* s using v -- change by MM */
ej = 2 * x * v;
v = v * v;
for (j = 1;; j++) { /* Taylor series */
ej *= v;
s1 = s + ej / ((j << 1) + 1);
if (s1 == s) /* last term was effectively 0 */ {
return (s1);
}
s = s1;
}
}
/* else: | x - np | is not too small */
return (x * Math.log(x / np) + np - x);
}
示例3: pnf
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double pnf(double x, double df1, double df2, double ncp, boolean lower_tail, boolean log_p) {
double y;
if (DoubleVector.isNaN(x) || DoubleVector.isNaN(df1) || DoubleVector.isNaN(df2) || DoubleVector.isNaN(ncp)) {
return x + df2 + df1 + ncp;
}
if (df1 <= 0. || df2 <= 0. || ncp < 0) {
return DoubleVector.NaN;
}
if (!DoubleVector.isFinite(ncp)) {
return DoubleVector.NaN;
}
if (!DoubleVector.isFinite(df1) && !DoubleVector.isFinite(df2)) { /* both +Inf */
return DoubleVector.NaN;
}
//R_P_bounds_01(x, 0., ML_POSINF);
if (x <= 0.0) {
return SignRank.R_DT_0(lower_tail, log_p);
}
if (x >= Double.POSITIVE_INFINITY) {
return SignRank.R_DT_1(lower_tail, log_p);
}
if (df2 > 1e8) /* avoid problems with +Inf and loss of accuracy */ {
return Distributions.pnchisq(x * df1, df1, ncp, lower_tail, log_p);
}
y = (df1 / df2) * x;
return Beta.pnbeta2(y / (1. + y), 1. / (1. + y), df1 / 2., df2 / 2.,
ncp, lower_tail, log_p);
}
示例4: qnf
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double qnf(double p, double df1, double df2, double ncp, boolean lower_tail, boolean log_p) {
double y;
if (DoubleVector.isNaN(p) || DoubleVector.isNaN(df1) || DoubleVector.isNaN(df2) || DoubleVector.isNaN(ncp)) {
return p + df1 + df2 + ncp;
}
if (df1 <= 0. || df2 <= 0. || ncp < 0) {
return DoubleVector.NaN;
}
if (!DoubleVector.isFinite(ncp)) {
return DoubleVector.NaN;
}
if (!DoubleVector.isFinite(df1) && !DoubleVector.isFinite(df2)) {
return DoubleVector.NaN;
}
//R_Q_P01_boundaries(p, 0, ML_POSINF);
if ((log_p && p > 0) || (!log_p && (p < 0 || p > 1))) {
return DoubleVector.NaN;
}
if (p == SignRank.R_DT_0(lower_tail, log_p)) {
return 0;
}
if (p == SignRank.R_DT_1(lower_tail, log_p)) {
return Double.POSITIVE_INFINITY;
}
//end of R_Q_P01_boundaries
if (df2 > 1e8) /* avoid problems with +Inf and loss of accuracy */ {
return Distributions.qnchisq(p, df1, ncp, lower_tail, log_p) / df1;
}
y = Beta.qnbeta(p, df1 / 2., df2 / 2., ncp, lower_tail, log_p);
return y / (1 - y) * (df2 / df1);
}
示例5: pnchisq
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double pnchisq(double x, double df, double ncp, boolean lower_tail, boolean log_p) {
double ans;
if (DoubleVector.isNaN(x) || DoubleVector.isNaN(df) || DoubleVector.isNaN(ncp)) {
return x + df + ncp;
}
if (!DoubleVector.isFinite(df) || !DoubleVector.isFinite(ncp)) {
return DoubleVector.NaN;
}
if (df < 0. || ncp < 0.) {
return DoubleVector.NaN;
}
ans = pnchisq_raw(x, df, ncp, 1e-12, 8 * SignRank.DBL_EPSILON, 1000000, lower_tail);
if (ncp >= 80) {
if (lower_tail) {
ans = Math.min(ans, 1.0); /* e.g., pchisq(555, 1.01, ncp = 80) */
} else { /* !lower_tail */
/* since we computed the other tail cancellation is likely */
if (ans < 1e-10) {
/*
* Is there any function for ML_ERROR?????
*/
//ML_ERROR(ME_PRECISION, "pnchisq");
}
ans = Math.max(ans, 0.0); /* Precaution PR#7099 */
}
}
if (!log_p) {
return ans;
}
/* if ans is near one, we can do better using the other tail */
if (ncp >= 80 || ans < 1 - 1e-8) {
return Math.log(ans);
}
ans = pnchisq_raw(x, df, ncp, 1e-12, 8 * SignRank.DBL_EPSILON, 1000000, !lower_tail);
return Math.log1p(-ans);
}
示例6: pgeom
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double pgeom(double x, double p, boolean lower_tail, boolean log_p) {
if (DoubleVector.isNaN(x) || DoubleVector.isNaN(p)) {
return x + p;
}
if (p <= 0 || p > 1) {
return DoubleVector.NaN;
}
if (x < 0.) {
return SignRank.R_DT_0(lower_tail, log_p);
}
if (!DoubleVector.isFinite(x)) {
return SignRank.R_DT_1(lower_tail, log_p);
}
x = Math.floor(x + 1e-7);
if (p == 1.) { /* we cannot assume IEEE */
x = lower_tail ? 1 : 0;
return log_p ? Math.log(x) : x;
}
x = Math.log1p(-p) * (x + 1);
if (log_p) {
return SignRank.R_DT_Clog(x, lower_tail, log_p);
} else {
return lower_tail ? -Math.expm1(x) : Math.exp(x);
}
}
示例7: dgeom
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double dgeom(double x, double p, boolean give_log) {
double prob;
if (DoubleVector.isNaN(x) || DoubleVector.isNaN(p)) {
return x + p;
}
if (p <= 0 || p > 1) {
return DoubleVector.NaN;
}
if (SignRank.R_D_nonint(x, true, give_log)) {
return SignRank.R_D__0(true, give_log);
}
if (x < 0 || !DoubleVector.isFinite(x) || p == 0) {
return SignRank.R_D__0(true, give_log);
}
x = SignRank.R_D_forceint(x);
/* prob = (1-p)^x, stable for small p */
prob = Binom.dbinom_raw(0., x, p, 1 - p, give_log);
return ((give_log) ? Math.log(p) + prob : p * prob);
}
示例8: pnbinom_mu
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double pnbinom_mu(double x, double size, double mu, boolean lower_tail, boolean log_p) {
if (DoubleVector.isNaN(x) || DoubleVector.isNaN(size) || DoubleVector.isNaN(mu)) {
return x + size + mu;
}
if (!DoubleVector.isFinite(size) || !DoubleVector.isFinite(mu)) {
return DoubleVector.NaN;
}
if (size <= 0 || mu < 0) {
return DoubleVector.NaN;
}
if (x < 0) {
SignRank.R_DT_0(lower_tail, log_p);
}
if (!DoubleVector.isFinite(x)) {
return SignRank.R_DT_1(lower_tail, log_p);
}
x = Math.floor(x + 1e-7);
/* return
* pbeta(pr, size, x + 1, lower_tail, log_p); pr = size/(size + mu), 1-pr = mu/(size+mu)
*
*= pbeta_raw(pr, size, x + 1, lower_tail, log_p)
* x. pin qin
*= bratio (pin, qin, x., 1-x., &w, &wc, &ierr, log_p), and return w or wc ..
*= bratio (size, x+1, pr, 1-pr, &w, &wc, &ierr, log_p) */
{
int[] ierr = new int[1];
double[] w = new double[1];
double[] wc = new double[1];
Utils.bratio(size, x + 1, size / (size + mu), mu / (size + mu), w, wc, ierr, log_p);
return lower_tail ? w[0] : wc[0];
}
}
示例9: dnf
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double dnf(double x, double df1, double df2, double ncp, boolean give_log) {
double y, z, f;
if (DoubleVector.isNaN(x) || DoubleVector.isNaN(df1) || DoubleVector.isNaN(df2) || DoubleVector.isNaN(ncp)) {
return x + df2 + df1 + ncp;
}
/* want to compare dnf(ncp=0) behavior with df() one, hence *NOT* :
* if (ncp == 0)
* return df(x, df1, df2, give_log); */
if (df1 <= 0. || df2 <= 0. || ncp < 0) {
return DoubleVector.NaN;
}
if (x < 0.) {
return (SignRank.R_D__0(true, give_log));
}
if (!DoubleVector.isFinite(ncp)) { /* ncp = +Inf -- FIXME?: in some cases, limit exists */
return DoubleVector.NaN;
}
/* This is not correct for df1 == 2, ncp > 0 - and seems unneeded:
* if (x == 0.) return(df1 > 2 ? R_D__0 : (df1 == 2 ? R_D__1 : ML_POSINF));
*/
if (!DoubleVector.isFinite(df1) && !DoubleVector.isFinite(df2)) { /* both +Inf */
/* PR: not sure about this (taken from ncp==0) -- FIXME ? */
if (x == 1.) {
return Double.POSITIVE_INFINITY;
}
/* else */ return SignRank.R_D__0(true, give_log);
}
if (!DoubleVector.isFinite(df2)) /* i.e. = +Inf */ {
return df1 * Distributions.dnchisq(x * df1, df1, ncp, give_log);
}
/* == dngamma(x, df1/2, 2./df1, ncp, give_log) -- but that does not exist */
if (df1 > 1e14 && ncp < 1e7) {
/* includes df1 == +Inf: code below is inaccurate there */
f = 1 + ncp / df1; /* assumes ncp << df1 [ignores 2*ncp^(1/2)/df1*x term] */
z = Distributions.dgamma(1. / x / f, df2 / 2, 2. / df2, give_log);
return give_log ? z - 2 * Math.log(x) - Math.log(f) : z / (x * x) / f;
}
y = (df1 / df2) * x;
z = Distributions.dnbeta(y / (1 + y), df1 / 2., df2 / 2., ncp, give_log);
return give_log
? z + Math.log(df1) - Math.log(df2) - 2 * Math.log1p(y)
: z * (df1 / df2) / (1 + y) / (1 + y);
}
示例10: rmultinom
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static void rmultinom(/*Session context,*/ int n, double[] prob, int K, int[] rN) /* `Return' vector rN[1:K] {K := length(prob)}
* where rN[j] ~ Bin(n, prob[j]) , sum_j rN[j] == n, sum_j prob[j] == 1,
*/ {
int k;
double pp;
double p_tot = 0.; /* LDOUBLE */
/* This calculation is sensitive to exact values, so we try to
ensure that the calculations are as accurate as possible
so different platforms are more likely to give the same
result. */
if (K < 1) {
//ML_ERROR(ME_DOMAIN, "rmultinom");
return;
}
if (n < 0) {
return;
}
if (K == IntVector.NA || K < 1) {
//ML_ERROR(ME_DOMAIN, "rmultinom");
return;
}
if (n == IntVector.NA || n < 0) {
return;
}
/* Note: prob[K] is only used here for checking sum_k prob[k] = 1 ;
* Could make loop one shorter and drop that check !
*/
for (k = 0; k < K; k++) {
pp = prob[k];
if (!DoubleVector.isFinite(pp) || pp < 0. || pp > 1.) {
//ML_ERR_ret_NAN(k);
return;
}
p_tot += pp;
rN[k] = 0;
}
if (Math.abs(p_tot - 1.) > 1e-7) //MATHLIB_ERROR(_("rbinom: probability sum should be 1, but is %g"), (double) p_tot);
{
if (n == 0) {
return;
}
}
if (K == 1 && p_tot == 0.) {
return; /* trivial border case: do as rbinom */
}
/* Generate the first K-1 obs. via binomials */
for (k = 0; k < K - 1; k++) { /* (p_tot, n) are for "remaining binomial" */
if (prob[k] != 0.0) {
pp = prob[k] / p_tot;
/* printf("[%d] %.17f\n", k+1, pp); */
rN[k] = ((pp < 1.) ? (int) Binom.rbinom(/*context,*/ (double) n, pp)
: /*>= 1; > 1 happens because of rounding */ n);
n -= rN[k];
} else {
rN[k] = 0;
}
if (n <= 0) /* we have all*/ {
return;
}
p_tot -= prob[k]; /* i.e. = sum(prob[(k+1):K]) */
}
rN[K - 1] = n;
return;
}
示例11: qnt
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double qnt(double p, double df, double ncp, boolean lower_tail, boolean log_p) {
final double accu = 1e-13;
final double Eps = 1e-11; /* must be > accu */
double ux, lx, nx, pp;
if (DoubleVector.isNaN(p) || DoubleVector.isNaN(df) || DoubleVector.isNaN(ncp)) {
return p + df + ncp;
}
if (!DoubleVector.isFinite(df)) {
return DoubleVector.NaN;
}
/* Was
* df = floor(df + 0.5);
* if (df < 1 || ncp < 0) ML_ERR_return_NAN;
*/
if (df <= 0.0) {
return (DoubleVector.NaN);
}
if (ncp == 0.0 && df >= 1.0) {
return Distributions.qt(p, df, lower_tail, log_p);
}
//R_Q_P01_boundaries(p, ML_NEGINF, ML_POSINF);
//#define R_Q_P01_boundaries(p, _LEFT_, _RIGHT_)
if (log_p) {
if (p > 0) {
return DoubleVector.NaN;
}
if (p == 0) /* upper bound*/ {
return lower_tail ? Double.POSITIVE_INFINITY : Double.NEGATIVE_INFINITY;
}
if (p == Double.NEGATIVE_INFINITY) {
return lower_tail ? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
}
} else { /* !log_p */
if (p < 0 || p > 1) {
return DoubleVector.NaN;
}
if (p == 0) {
return lower_tail ? Double.NEGATIVE_INFINITY : Double.POSITIVE_INFINITY;
}
if (p == 1) {
return lower_tail ? Double.POSITIVE_INFINITY : Double.NEGATIVE_INFINITY;
}
}
p = Normal.R_DT_qIv(p, log_p ? 1.0 : 0.0, lower_tail ? 1.0 : 0.0);
/* Invert pnt(.) :
* 1. finding an upper and lower bound */
if (p > 1 - SignRank.DBL_EPSILON) {
return Double.POSITIVE_INFINITY;
}
pp = Math.min(1 - SignRank.DBL_EPSILON, p * (1 + Eps));
for (ux = Math.max(1., ncp);
ux < Double.MAX_VALUE && pnt(ux, df, ncp, true, false) < pp;
ux *= 2);
pp = p * (1 - Eps);
for (lx = Math.min(-1., -ncp);
lx > -Double.MAX_VALUE && pnt(lx, df, ncp, true, false) > pp;
lx *= 2);
/* 2. interval (lx,ux) halving : */
do {
nx = 0.5 * (lx + ux);
if (pnt(nx, df, ncp, true, false) > p) {
ux = nx;
} else {
lx = nx;
}
} while ((ux - lx) / Math.abs(nx) > accu);
return 0.5 * (lx + ux);
}
示例12: dnt
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double dnt(double x, double df, double ncp, boolean give_log) {
double u;
if (DoubleVector.isNaN(x) || DoubleVector.isNaN(df)) {
return x + df;
}
/* If non-positive df then error */
if (df <= 0.0) {
return DoubleVector.NaN;
}
if (ncp == 0.0) {
return Distributions.dt(x, df, give_log);
}
/* If x is infinite then return 0 */
if (!DoubleVector.isFinite(x)) {
return SignRank.R_D__0(true, give_log);
}
/* If infinite df then the density is identical to a
normal distribution with mean = ncp. However, the formula
loses a lot of accuracy around df=1e9
*/
if (!DoubleVector.isFinite(df) || df > 1e8) {
return Distributions.dnorm(x, ncp, 1., give_log);
}
/* Do calculations on log scale to stabilize */
/* Consider two cases: x ~= 0 or not */
if (Math.abs(x) > Math.sqrt(df * SignRank.DBL_EPSILON)) {
u = Math.log(df) - Math.log(Math.abs(x))
+ Math.log(Math.abs(Distributions.pnt(x * Math.sqrt((df + 2) / df), df + 2, ncp, true, false)
- Distributions.pnt(x, df, ncp, true, false)));
/* FIXME: the above still suffers from cancellation (but not horribly) */
} else { /* x ~= 0 : -> same value as for x = 0 */
u = org.apache.commons.math.special.Gamma.logGamma((df + 1) / 2) - org.apache.commons.math.special.Gamma.logGamma(df / 2)
- .5 * (Math.log(Math.PI) + Math.log(df) + ncp * ncp);
}
return (give_log ? u : Math.exp(u));
}
示例13: dnbeta
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double dnbeta(double x, double a, double b, double ncp, boolean give_log) {
final double eps = 1.e-15;
int kMax;
double k, ncp2, dx2, d, D;
double sum, term, p_k, q; /* They were LDOUBLE */
if (DoubleVector.isNaN(x) || DoubleVector.isNaN(a) || DoubleVector.isNaN(b) || DoubleVector.isNaN(ncp)) {
return x + a + b + ncp;
}
if (ncp < 0 || a <= 0 || b <= 0) {
return DoubleVector.NaN;
}
if (!DoubleVector.isFinite(a) || !DoubleVector.isFinite(b) || !DoubleVector.isFinite(ncp)) {
return DoubleVector.NaN;
}
if (x < 0 || x > 1) {
return (SignRank.R_D__0(true, give_log));
}
if (ncp == 0) {
return Distributions.dbeta(x, a, b, give_log);
}
/* New algorithm, starting with *largest* term : */
ncp2 = 0.5 * ncp;
dx2 = ncp2 * x;
d = (dx2 - a - 1) / 2;
D = d * d + dx2 * (a + b) - a;
if (D <= 0) {
kMax = 0;
} else {
D = Math.ceil(d + Math.sqrt(D));
kMax = (D > 0) ? (int) D : 0;
}
/* The starting "middle term" --- first look at it's log scale: */
term = Distributions.dbeta(x, a + kMax, b, /* log = */ true);
p_k = Poisson.dpois_raw(kMax, ncp2, true);
if (x == 0. || !DoubleVector.isFinite(term) || !DoubleVector.isFinite(p_k)) /* if term = +Inf */ {
return SignRank.R_D_exp(p_k + term, true, give_log);
}
/* Now if s_k := p_k * t_k {here = exp(p_k + term)} would underflow,
* we should rather scale everything and re-scale at the end:*/
p_k += term; /* = log(p_k) + log(t_k) == log(s_k) -- used at end to rescale */
/* mid = 1 = the rescaled value, instead of mid = exp(p_k); */
/* Now sum from the inside out */
sum = term = 1. /* = mid term */;
/* middle to the left */
k = kMax;
while (k > 0 && term > sum * eps) {
k--;
q = /* 1 / r_k = */ (k + 1) * (k + a) / (k + a + b) / dx2;
term *= q;
sum += term;
}
/* middle to the right */
term = 1.;
k = kMax;
do {
q = /* r_{old k} = */ dx2 * (k + a + b) / (k + a) / (k + 1);
k++;
term *= q;
sum += term;
} while (term > sum * eps);
return SignRank.R_D_exp(p_k + Math.log(sum), true, give_log);
}
示例14: qnbeta
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double qnbeta(double p, double a, double b, double ncp,
/*@Recycle(false)*/ boolean lower_tail, /*@Recycle(false)*/ boolean log_p) {
final double accu = 1e-15;
final double Eps = 1e-14; /* must be > accu */
double ux, lx, nx, pp;
if (DoubleVector.isNaN(p) || DoubleVector.isNaN(a) || DoubleVector.isNaN(b) || DoubleVector.isNaN(ncp)) {
return p + a + b + ncp;
}
if (!DoubleVector.isFinite(a)) {
return DoubleVector.NaN;
}
if (ncp < 0. || a <= 0. || b <= 0.) {
return DoubleVector.NaN;
}
//R_Q_P01_boundaries(p, 0, 1);
if ((log_p && p > 0) || (!log_p && (p < 0 || p > 1))) {
return DoubleVector.NaN;
}
if (p == SignRank.R_DT_0(lower_tail, log_p)) {
return 0.0;
}
if (p == SignRank.R_DT_1(lower_tail, log_p)) {
return 1.0;
}
//end of R_Q_P01_boundaries
p = Normal.R_DT_qIv(p, log_p ? 1.0 : 0.0, lower_tail ? 1.0 : 0.0);
/* Invert pnbeta(.) :
* 1. finding an upper and lower bound */
if (p > 1 - SignRank.DBL_EPSILON) {
return 1.0;
}
pp = Math.min(1 - SignRank.DBL_EPSILON, p * (1 + Eps));
for (ux = 0.5;
ux < 1 - SignRank.DBL_EPSILON && pnbeta(ux, a, b, ncp, true, false) < pp;
ux = 0.5 * (1 + ux));
pp = p * (1 - Eps);
for (lx = 0.5;
lx > Double.MIN_VALUE && pnbeta(lx, a, b, ncp, true, false) > pp;
lx *= 0.5);
/* 2. interval (lx,ux) halving : */
do {
nx = 0.5 * (lx + ux);
if (pnbeta(nx, a, b, ncp, true, false) > p) {
ux = nx;
} else {
lx = nx;
}
} while ((ux - lx) / nx > accu);
return 0.5 * (ux + lx);
}
示例15: dnbinom_mu
import org.renjin.sexp.DoubleVector; //导入方法依赖的package包/类
public static double dnbinom_mu(double x, double size, double mu, boolean give_log) {
/* originally, just set prob := size / (size + mu) and called dbinom_raw(),
* but that suffers from cancellation when mu << size */
double ans, p;
if (DoubleVector.isNaN(x) || DoubleVector.isNaN(size) || DoubleVector.isNaN(mu)) {
return x + size + mu;
}
if (mu < 0 || size < 0) {
return DoubleVector.NaN;
}
//R_D_nonint_check(x);
if (SignRank.R_D_nonint(x, true, give_log)) {
//MATHLIB_WARNING("non-integer x = %f", x);
//How to warn??
return SignRank.R_D__0(true, give_log);
}
if (x < 0 || !DoubleVector.isFinite(x)) {
return SignRank.R_D__0(true, give_log);
}
x = SignRank.R_D_forceint(x);
if (x == 0)/* be accerate, both for n << mu, and n >> mu :*/ {
return SignRank.R_D_exp(size * (size < mu ? Math.log(size / (size + mu)) : Math.log1p(-mu / (size + mu))), true, give_log);
}
if (x < 1e-10 * size) { /* don't use dbinom_raw() but MM's formula: */
/* FIXME --- 1e-8 shows problem; rather use algdiv() from ./toms708.c */
return SignRank.R_D_exp(x * Math.log(size * mu / (size + mu)) - mu - org.apache.commons.math.special.Gamma.logGamma(x + 1)
+ Math.log1p(x * (x - 1) / (2 * size)), true, give_log);
}
/* else: no unnecessary cancellation inside dbinom_raw, when
* x_ = size and n_ = x+size are so close that n_ - x_ loses accuracy
*/
ans = dbinom_raw(size, x + size, size / (size + mu), mu / (size + mu), give_log);
p = ((double) size) / (size + x);
return ((give_log) ? Math.log(p) + ans : p * ans);
}