Skip to content
Snippets Groups Projects
Commit b52dac65 authored by stepan's avatar stepan
Browse files

Implement {x}lnorm functions and df

parent f83c42a2
No related branches found
No related tags found
No related merge requests found
/*
* This material is distributed under the GNU General Public License
* Version 2. You may review the terms of this license at
* http://www.gnu.org/licenses/gpl-2.0.html
*
* Copyright (c) 2000--2014, The R Core Team
* Copyright (c) 2016, Oracle and/or its affiliates
*
* All rights reserved.
*/
// Acknowledgement from GnuR header:
// Author: Catherine Loader, catherine@research.bell-labs.com, October 23, 2000.
package com.oracle.truffle.r.library.stats;
import static com.oracle.truffle.r.library.stats.Dbinom.dbinomRaw;
import static com.oracle.truffle.r.library.stats.GammaFunctions.dgamma;
import com.oracle.truffle.r.library.stats.StatsFunctions.Function3_1;
public final class Df implements Function3_1 {
@Override
public double evaluate(double x, double m, double n, boolean giveLog) {
double p;
double q;
double f;
double dens;
if (Double.isNaN(x) || Double.isNaN(m) || Double.isNaN(n)) {
return x + m + n;
}
if (m <= 0 || n <= 0) {
return RMath.mlError();
}
if (x < 0.) {
return DPQ.rd0(giveLog);
}
if (x == 0.) {
return m > 2 ? DPQ.rd0(giveLog) : (m == 2 ? DPQ.rd1(giveLog) : Double.POSITIVE_INFINITY);
}
if (!Double.isFinite(m) && !Double.isFinite(n)) { /* both +Inf */
if (x == 1.) {
return Double.POSITIVE_INFINITY;
} else {
return DPQ.rd0(giveLog);
}
}
if (!Double.isFinite(n)) {
/* must be +Inf by now */
return dgamma(x, m / 2, 2. / m, giveLog);
}
if (m > 1e14) { /* includes +Inf: code below is inaccurate there */
dens = dgamma(1. / x, n / 2, 2. / n, giveLog);
return giveLog ? dens - 2 * Math.log(x) : dens / (x * x);
}
f = 1. / (n + x * m);
q = n * f;
p = x * m * f;
if (m >= 2) {
f = m * q / 2;
dens = dbinomRaw((m - 2) / 2, (m + n - 2) / 2, p, q, giveLog);
} else {
f = m * m * q / (2 * p * (m + n));
dens = dbinomRaw(m / 2, (m + n) / 2, p, q, giveLog);
}
return (giveLog ? Math.log(f) + dens : f * dens);
}
}
/*
* This material is distributed under the GNU General Public License
* Version 2. You may review the terms of this license at
* http://www.gnu.org/licenses/gpl-2.0.html
*
* Copyright (C) 1998 Ross Ihaka
* Copyright (c) 2000--2014, The R Core Team
* Copyright (c) 2005, The R Foundation
* Copyright (c) 2016, Oracle and/or its affiliates
*
* All rights reserved.
*/
package com.oracle.truffle.r.library.stats;
import com.oracle.truffle.r.library.stats.DPQ.EarlyReturn;
import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandFunction2_Double;
import com.oracle.truffle.r.library.stats.RandGenerationFunctions.RandomNumberProvider;
import com.oracle.truffle.r.library.stats.StatsFunctions.Function3_1;
import com.oracle.truffle.r.library.stats.StatsFunctions.Function3_2;
public final class LogNormal {
private LogNormal() {
// only static members
}
public static final class RLNorm implements RandFunction2_Double {
private final Rnorm rnorm = new Rnorm();
@Override
public double evaluate(double meanlog, double sdlog, RandomNumberProvider rand) {
if (Double.isNaN(meanlog) || !Double.isFinite(sdlog) || sdlog < 0.) {
return RMath.mlError();
}
return Math.exp(rnorm.evaluate(meanlog, sdlog, rand));
}
}
public static final class DLNorm implements Function3_1 {
@Override
public double evaluate(double x, double meanlog, double sdlog, boolean giveLog) {
if (Double.isNaN(x) || Double.isNaN(meanlog) || Double.isNaN(sdlog)) {
return x + meanlog + sdlog;
}
if (sdlog <= 0) {
if (sdlog < 0) {
return RMath.mlError();
}
// sdlog == 0 :
return (Math.log(x) == meanlog) ? Double.POSITIVE_INFINITY : DPQ.rd0(giveLog);
}
if (x <= 0) {
return DPQ.rd0(giveLog);
}
double y = (Math.log(x) - meanlog) / sdlog;
return (giveLog ? -(MathConstants.M_LN_SQRT_2PI + 0.5 * y * y + Math.log(x * sdlog)) : MathConstants.M_1_SQRT_2PI * Math.exp(-0.5 * y * y) / (x * sdlog));
/* M_1_SQRT_2PI = 1 / Math.sqrt(2 * pi) */
}
}
public static final class QLNorm implements Function3_2 {
private final Qnorm qnorm = new Qnorm();
@Override
public double evaluate(double p, double meanlog, double sdlog, boolean lowerTail, boolean logP) {
if (Double.isNaN(p) || Double.isNaN(meanlog) || Double.isNaN(sdlog)) {
return p + meanlog + sdlog;
}
try {
DPQ.rqp01boundaries(p, 0, Double.POSITIVE_INFINITY, lowerTail, logP);
} catch (EarlyReturn e) {
return e.result;
}
return Math.exp(qnorm.evaluate(p, meanlog, sdlog, lowerTail, logP));
}
}
public static final class PLNorm implements Function3_2 {
private final Pnorm pnorm = new Pnorm();
@Override
public double evaluate(double x, double meanlog, double sdlog, boolean lowerTail, boolean logP) {
if (Double.isNaN(x) || Double.isNaN(meanlog) || Double.isNaN(sdlog)) {
return x + meanlog + sdlog;
}
if (sdlog < 0) {
return RMath.mlError();
}
if (x > 0) {
return pnorm.evaluate(Math.log(x), meanlog, sdlog, lowerTail, logP);
}
return DPQ.rdt0(lowerTail, logP);
}
}
}
......@@ -41,11 +41,11 @@ public final class Pnorm implements StatsFunctions.Function3_2 {
return Double.NaN;
}
/* sigma = 0 : */
return (x < mu) ? DPQ.rd0(logP) : DPQ.rd1(logP);
return (x < mu) ? DPQ.rdt0(lowerTail, logP) : DPQ.rdt1(lowerTail, logP);
}
double p = (x - mu) / sigma;
if (!Double.isFinite(p)) {
return (x < mu) ? DPQ.rd0(logP) : DPQ.rd1(logP);
return (x < mu) ? DPQ.rdt0(lowerTail, logP) : DPQ.rdt1(lowerTail, logP);
}
PnormBoth pnormBoth = new PnormBoth(p);
......
......@@ -52,6 +52,7 @@ import com.oracle.truffle.r.library.stats.CutreeNodeGen;
import com.oracle.truffle.r.library.stats.DBeta;
import com.oracle.truffle.r.library.stats.DPois;
import com.oracle.truffle.r.library.stats.Dbinom;
import com.oracle.truffle.r.library.stats.Df;
import com.oracle.truffle.r.library.stats.DoubleCentreNodeGen;
import com.oracle.truffle.r.library.stats.Dt;
import com.oracle.truffle.r.library.stats.Exp.DExp;
......@@ -63,6 +64,10 @@ import com.oracle.truffle.r.library.stats.GammaFunctions.QgammaFunc;
import com.oracle.truffle.r.library.stats.Geom;
import com.oracle.truffle.r.library.stats.Geom.DGeom;
import com.oracle.truffle.r.library.stats.Geom.RGeom;
import com.oracle.truffle.r.library.stats.LogNormal;
import com.oracle.truffle.r.library.stats.LogNormal.DLNorm;
import com.oracle.truffle.r.library.stats.LogNormal.PLNorm;
import com.oracle.truffle.r.library.stats.LogNormal.QLNorm;
import com.oracle.truffle.r.library.stats.Pbeta;
import com.oracle.truffle.r.library.stats.Pbinom;
import com.oracle.truffle.r.library.stats.Pf;
......@@ -299,6 +304,8 @@ public class CallAndExternalFunctions {
return StatsFunctionsFactory.Function3_2NodeGen.create(new Cauchy.QCauchy());
case "pf":
return StatsFunctionsFactory.Function3_2NodeGen.create(new Pf());
case "df":
return StatsFunctionsFactory.Function3_1NodeGen.create(new Df());
case "dgamma":
return StatsFunctionsFactory.Function3_1NodeGen.create(new DGamma());
case "dchisq":
......@@ -321,6 +328,14 @@ public class CallAndExternalFunctions {
return StatsFunctionsFactory.Function3_1NodeGen.create(new DBeta());
case "dt":
return StatsFunctionsFactory.Function2_1NodeGen.create(new Dt());
case "rlnorm":
return RandGenerationFunctionsFactory.Function2_DoubleNodeGen.create(new LogNormal.RLNorm());
case "dlnorm":
return StatsFunctionsFactory.Function3_1NodeGen.create(new DLNorm());
case "qlnorm":
return StatsFunctionsFactory.Function3_2NodeGen.create(new QLNorm());
case "plnorm":
return StatsFunctionsFactory.Function3_2NodeGen.create(new PLNorm());
case "rmultinom":
return RMultinomNodeGen.create();
case "Approx":
......
......@@ -31,7 +31,7 @@ import com.oracle.truffle.r.test.TestBase;
* tests for its specific corner cases if those are not covered here.
*/
public class TestRandGenerationFunctions extends TestBase {
private static final String[] FUNCTION2_NAMES = {"rnorm", "runif", "rgamma", "rbeta", "rcauchy", "rf", "rlogis", "rweibull", "rchisq", "rwilcox"};
private static final String[] FUNCTION2_NAMES = {"rnorm", "runif", "rgamma", "rbeta", "rcauchy", "rf", "rlogis", "rweibull", "rchisq", "rwilcox", "rlnorm"};
private static final String[] FUNCTION2_PARAMS = {
"10, 10, 10",
"20, c(-1, 0, 0.2, 2:5), c(-1, 0, 0.1, 0.9, 3)",
......
......@@ -30,7 +30,7 @@ import com.oracle.truffle.r.test.TestBase;
* Common tests for functions implemented using {@code StatsFunctions} infrastructure.
*/
public class TestStatFunctions extends TestBase {
private static final String[] FUNCTION3_1_NAMES = {"dgamma", "dbeta", "dcauchy"};
private static final String[] FUNCTION3_1_NAMES = {"dgamma", "dbeta", "dcauchy", "dlnorm"};
private static final String[] FUNCTION3_1_PARAMS = {
"10, 10, 10, log=TRUE",
"3, 3, 3, log=FALSE",
......@@ -80,7 +80,7 @@ public class TestStatFunctions extends TestBase {
assertEval(Output.IgnoreWhitespace, template("set.seed(1); %0(%1)", FUNCTION2_2_NAMES, new String[]{"rep(c(1, 0, 0.1), 5), c(NA, 0, NaN, 1/0, -1/0)"}));
}
private static final String[] FUNCTION3_2_NAMES = {"pbeta", "pcauchy", "qcauchy"};
private static final String[] FUNCTION3_2_NAMES = {"pbeta", "pcauchy", "qcauchy", "qlnorm", "plnorm"};
private static final String[] FUNCTION3_2_PARAMS = {
"0, 10, 10",
"c(-1, 0, 0.2, 2), c(-1, 0, 0.1, 0.9, 3), rep(c(-1, 0, 1, 0.1, -0.1, 0.0001), 20)",
......
......@@ -61,6 +61,7 @@ com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/SplineFuncti
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/StatsFunctions.java,gnu_r_gentleman_ihaka.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RandGenerationFunctions.java,gnu_r_gentleman_ihaka.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RMath.java,gnu_r_ihaka.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/LogNormal.java,gnu_r_ihaka.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/MathConstants.java,gnu_r_ihaka.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/TOMS708.java,gnu_r.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/SNorm.java,gnu_r_ihaka_core.copyright
......@@ -73,6 +74,7 @@ com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RChisq.java,
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Exp.java,gnu_r_ihaka_core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Geom.java,gnu_r_ihaka_core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Dt.java,gnu_r.core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/Df.java,gnu_r.core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DBeta.java,gnu_r.core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/DPois.java,gnu_r.core.copyright
com.oracle.truffle.r.library/src/com/oracle/truffle/r/library/stats/RNchisq.java,gnu_r.copyright
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment