List of usage examples for java.lang Math exp
@HotSpotIntrinsicCandidate public static double exp(double a)
From source file:com.datumbox.framework.statistics.distributions.ContinuousDistributions.java
/** * Internal function used by GammaCdf/*w w w .j av a2 s. c om*/ * * @param x * @param a * @return * @throws IllegalArgumentException */ protected static double GammaCdf(double x, double a) throws IllegalArgumentException { if (x < 0) { throw new IllegalArgumentException(); } double GI = 0; if (a > 200) { double z = (x - a) / Math.sqrt(a); double y = GaussCdf(z); double b1 = 2 / Math.sqrt(a); double phiz = 0.39894228 * Math.exp(-z * z / 2); double w = y - b1 * (z * z - 1) * phiz / 6; //Edgeworth1 double b2 = 6 / a; int zXor4 = ((int) z) ^ 4; double u = 3 * b2 * (z * z - 3) + b1 * b1 * (zXor4 - 10 * z * z + 15); GI = w - phiz * z * u / 72; //Edgeworth2 } else if (x < a + 1) { GI = Gser(x, a); } else { GI = Gcf(x, a); } return GI; }
From source file:edu.byu.nlp.util.Matrices.java
public static double logSumRowInColumnMajorMatrix(double[] mat, int numRows, int row) { double sumExponentiatedDiffs = 0.0; int argMax = rowArgMaxInColumnMajorMatrix(mat, numRows, row); double max = mat[argMax]; for (int i = row; i < argMax; i += numRows) { double diff = mat[i] - max; if (diff >= Math2.LOGSUM_THRESHOLD) { sumExponentiatedDiffs += Math.exp(diff); }// w w w . jav a 2s.c o m } for (int i = argMax + numRows; i < mat.length; i += numRows) { double diff = mat[i] - max; if (diff >= Math2.LOGSUM_THRESHOLD) { sumExponentiatedDiffs += Math.exp(diff); } } return max + Math.log(1.0 + sumExponentiatedDiffs); }
From source file:com.opengamma.analytics.financial.model.finitedifference.ForwardPDETest.java
@Test public void flatVolTest() { final double spot = 10.0; final double r = 0.1; final double y = 0.0; final double vol = 0.3; final double expiry = 2.0; final double mult = Math.exp(Math.sqrt(expiry) * vol * 6.0); final ConvectionDiffusionPDESolver solver = new ThetaMethodFiniteDifference(0.5, true); // ConvectionDiffusionPDESolver solver = new RichardsonExtrapolationFiniteDifference(base); final int tNodes = 51; final int xNodes = 101; final ConvectionDiffusionPDE1DStandardCoefficients pde = PDE_DATA_PROVIDER.getForwardBlackSholes(r, y, vol); final Function1D<Double, Double> intCon = INITIAL_COND_PROVIDER.getForwardCallPut(spot, true); //only true if y =0 final BoundaryCondition lower = new DirichletBoundaryCondition(spot, 0); final BoundaryCondition upper = new DirichletBoundaryCondition(0, mult * spot); final MeshingFunction timeMesh = new ExponentialMeshing(0, expiry, tNodes, 3.0); //final MeshingFunction spaceMesh = new ExponentialMeshing(0, 5 * spot, xNodes, 0.0); final MeshingFunction spaceMesh = new HyperbolicMeshing(0, mult * spot, spot, xNodes, 0.05); final PDEGrid1D grid = new PDEGrid1D(timeMesh, spaceMesh); // final PDE1DDataBundle<ConvectionDiffusionPDE1DCoefficients> db = new PDE1DDataBundle<ConvectionDiffusionPDE1DCoefficients>(pde, intCon, lower, upper, grid); // final PDEFullResults1D res = (PDEFullResults1D) solver.solve(db); final double[][] price = new double[tNodes][xNodes]; for (int i = 0; i < tNodes; i++) { final double t = grid.getTimeNode(i); final double df = Math.exp(-r * t); final double fwd = spot * Math.exp((r - y) * t); for (int j = 0; j < xNodes; j++) { price[i][j] = df * BlackFormulaRepository.price(fwd, grid.getSpaceNode(j), t, vol, true); }/*w w w .ja v a 2 s. com*/ } // double k = grid.getSpaceNode(1); // double pdePrice = res.getFunctionValue(1, 1); // double analPrice = price[1][1]; // double t = grid.getTimeNode(1); // double fPrice = pdePrice * Math.exp(r * t); // double fwd = spot * Math.exp(r * t); // double iv = BlackFormulaRepository.impliedVolatility(fPrice, fwd, k, t, true); // // System.out.println("debug " + grid.getSpaceNode(1) + "\t" + price[1][1] + "\t" + res.getFunctionValue(1, 1) + "\t" + iv); // PDEUtilityTools.printSurface("PDE res", res); // PDEUtilityTools.printSurface("call price", price, grid.getSpaceNodes(), grid.getTimeNodes(), System.out); // // Map<DoublesPair, Double> iv = PDEUtilityTools.priceToImpliedVol(new ForwardCurve(spot, r), new YieldCurve("", ConstantDoublesCurve.from(r)), res, 0, 2.0, 0.0, mult * spot); // GridInterpolator2D gridIn = new GridInterpolator2D(Interpolator1DFactory.DOUBLE_QUADRATIC_INSTANCE, Interpolator1DFactory.DOUBLE_QUADRATIC_INSTANCE); // Map<Double, Interpolator1DDataBundle> idb = gridIn.getDataBundle(iv); // PDEUtilityTools.printSurface("iv", idb, 0.1, 2.0, 0.7, 3 * spot, 100, 100); }
From source file:StatFunctions.java
public double pchisq(double q, double df) { // Posten, H. (1989) American Statistician 43 p. 261-265 double df2 = df * .5; double q2 = q * .5; int n = 5, k; double tk, CFL, CFU, prob; if (q <= 0 || df <= 0) throw new IllegalArgumentException("Illegal argument " + q + " or " + df + " for qnorm(p)."); if (q < df) { tk = q2 * (1 - n - df2) / (df2 + 2 * n - 1 + n * q2 / (df2 + 2 * n)); for (k = n - 1; k > 1; k--) tk = q2 * (1 - k - df2) / (df2 + 2 * k - 1 + k * q2 / (df2 + 2 * k + tk)); CFL = 1 - q2 / (df2 + 1 + q2 / (df2 + 2 + tk)); prob = Math.exp(df2 * Math.log(q2) - q2 - logGamma(df2 + 1) - Math.log(CFL)); } else {// w ww. j a v a 2 s .c om tk = (n - df2) / (q2 + n); for (k = n - 1; k > 1; k--) tk = (k - df2) / (q2 + k / (1 + tk)); CFU = 1 + (1 - df2) / (q2 + 1 / (1 + tk)); prob = 1 - Math.exp((df2 - 1) * Math.log(q2) - q2 - logGamma(df2) - Math.log(CFU)); } return prob; }
From source file:com.opengamma.analytics.financial.model.option.pricing.tree.NormalBinomialTreeBuilderTest.java
@Test public void testCEV() { final GeneralNormalOptionDataBundle data = new GeneralNormalOptionDataBundle(YIELD_CURVE, DRIFTLESS, new VolatilitySurface(FunctionalDoublesSurface.from(CEV_LOCAL_VOL)), FORWARD, DATE); final RecombiningBinomialTree<BinomialTreeNode<Double>> assetPriceTree = BUILDER.buildAssetTree(T, data, 200);// www.j a v a 2 s . c o m RecombiningBinomialTree<BinomialTreeNode<Double>> optionPriceTree = BUILDER.buildOptionPriceTree(OPTION, data, assetPriceTree); for (int i = 0; i < 10; i++) { final double m = -1.5 + 3.0 * i / 10.0; final double strike = FORWARD * Math.exp(ATM_VOL * Math.sqrt(T) * m); final OptionDefinition option = new EuropeanVanillaOptionDefinition(strike, OPTION.getExpiry(), OPTION.isCall()); optionPriceTree = BUILDER.buildOptionPriceTree(option, data, assetPriceTree); final EuropeanVanillaOption o = new EuropeanVanillaOption(strike, T, true); final CEVFunctionData cfd = new CEVFunctionData(FORWARD, YIELD_CURVE.getDiscountFactor(T), SIGMA_BETA, BETA); final double cevPrice = CEV_PRICE.getPriceFunction(o).evaluate(cfd); final double cevVol = BLACK_IMPLIED_VOL.getImpliedVolatility( new BlackFunctionData(FORWARD, YIELD_CURVE.getDiscountFactor(T), SIGMA_BETA), o, cevPrice); final double impVol = BLACK_IMPLIED_VOL.getImpliedVolatility( new BlackFunctionData(FORWARD, YIELD_CURVE.getDiscountFactor(T), SIGMA_BETA), o, optionPriceTree.getNode(0, 0).getValue()); //final double cevPrice = CEVFormula.optionPrice(FORWARD, strike, BETA, df, SIGMA_BETA, T, true); //final double cevVol = BlackImpliedVolFormula.impliedVol(cevPrice, FORWARD, strike, df, T, true); //final double impVol = BlackImpliedVolFormula.impliedVol(optionPriceTree.getNode(0, 0).getValue(), FORWARD, strike, df, T, true); // System.out.println(strike + "\t" + cevVol + "\t" + impVol); assertEquals(cevVol, impVol, 1e-3); } }
From source file:com.opengamma.analytics.financial.interestrate.swaption.provider.SwaptionPhysicalFixedIborLMMDDMethod.java
/** * Computes the present value of the Physical delivery swaption. * @param swaption The swaption.//from w w w .ja v a2s . c om * @param lmmData The LMM and multi-curves provider. * @return The present value. */ public MultipleCurrencyAmount presentValue(final SwaptionPhysicalFixedIbor swaption, final LiborMarketModelDisplacedDiffusionProviderInterface lmmData) { ArgumentChecker.notNull(swaption, "Swaption"); ArgumentChecker.notNull(lmmData, "LMM provider"); final Currency ccy = swaption.getCurrency(); final MulticurveProviderInterface multicurves = lmmData.getMulticurveProvider(); final LiborMarketModelDisplacedDiffusionParameters parameters = lmmData.getLMMParameters(); // 1. Swaption CFE preparation final AnnuityPaymentFixed cfe = swaption.getUnderlyingSwap().accept(CFEC, multicurves); final int nbCFInit = cfe.getNumberOfPayments(); final double multFact = Math.signum(cfe.getNthPayment(0).getAmount()); final boolean isCall = (cfe.getNthPayment(0).getAmount() < 0); final double[] cftInit = new double[nbCFInit]; final double[] cfaInit = new double[nbCFInit]; for (int loopcf = 0; loopcf < nbCFInit; loopcf++) { cftInit[loopcf] = cfe.getNthPayment(loopcf).getPaymentTime(); cfaInit[loopcf] = cfe.getNthPayment(loopcf).getAmount() * -multFact; } final double timeToExpiry = swaption.getTimeToExpiry(); // 2. Model data final int nbFactor = parameters.getNbFactor(); final double[][] volLMM = parameters.getVolatility(); final double[] timeLMM = parameters.getIborTime(); // 3. Link cfe dates to lmm final int[] indCFDate = new int[nbCFInit]; int indStart = nbCFInit - 1; int indEnd = 0; for (int loopcf = 0; loopcf < nbCFInit; loopcf++) { indCFDate[loopcf] = Arrays.binarySearch(timeLMM, cftInit[loopcf]); if (indCFDate[loopcf] < 0) { if (timeLMM[-indCFDate[loopcf] - 1] - cftInit[loopcf] < TIME_TOLERANCE) { indCFDate[loopcf] = -indCFDate[loopcf] - 1; } else { if (cftInit[loopcf] - timeLMM[-indCFDate[loopcf] - 2] < TIME_TOLERANCE) { indCFDate[loopcf] = -indCFDate[loopcf] - 2; } else { Validate.isTrue(true, "Instrument time incompatible with LMM parametrisation"); } } } if (indCFDate[loopcf] < indStart) { indStart = indCFDate[loopcf]; } if (indCFDate[loopcf] > indEnd) { indEnd = indCFDate[loopcf]; } } final int nbCF = indEnd - indStart + 1; final double[] cfa = new double[nbCF]; for (int loopcf = 0; loopcf < nbCFInit; loopcf++) { cfa[indCFDate[loopcf] - indStart] = cfaInit[loopcf]; } final double[] cft = new double[nbCF]; System.arraycopy(timeLMM, indStart, cft, 0, nbCF); final double[] dfLMM = new double[nbCF]; for (int loopcf = 0; loopcf < nbCF; loopcf++) { dfLMM[loopcf] = multicurves.getDiscountFactor(ccy, cft[loopcf]); } final double[][] gammaLMM = new double[nbCF - 1][nbFactor]; final double[] deltaLMM = new double[nbCF - 1]; System.arraycopy(parameters.getAccrualFactor(), indStart, deltaLMM, 0, nbCF - 1); final double[] aLMM = new double[nbCF - 1]; System.arraycopy(parameters.getDisplacement(), indStart, aLMM, 0, nbCF - 1); final double[] liborLMM = new double[nbCF - 1]; final double amr = parameters.getMeanReversion(); for (int loopcf = 0; loopcf < nbCF - 1; loopcf++) { gammaLMM[loopcf] = volLMM[indStart + loopcf]; liborLMM[loopcf] = (dfLMM[loopcf] / dfLMM[loopcf + 1] - 1.0d) / deltaLMM[loopcf]; } // TODO: 4. cfe modification (for roller coasters) final double[] cfaMod = new double[nbCF + 1]; final double cfaMod0 = cfa[0]; cfaMod[0] = cfaMod0; // modified strike cfaMod[1] = 0.0; System.arraycopy(cfa, 1, cfaMod, 2, nbCF - 1); // 5. Pricing algorithm final double[] p0 = new double[nbCF]; final double[] dP = new double[nbCF]; double b0 = 0; for (int loopcf = 0; loopcf < nbCF; loopcf++) { p0[loopcf] = dfLMM[loopcf] / dfLMM[0]; dP[loopcf] = cfaMod[loopcf + 1] * p0[loopcf]; b0 += dP[loopcf]; } final double bK = -cfaMod0; final double bM = (b0 + bK) / 2.0d; final double meanReversionImpact = Math.abs(amr) < 1.0E-6 ? timeToExpiry : (Math.exp(2.0d * amr * timeToExpiry) - 1.0d) / (2.0d * amr); // To handle 0 mean reversion. final double[] rate0Ratio = new double[nbCF - 1]; final double[][] mu0 = new double[nbCF - 1][nbFactor]; for (int loopcf = 0; loopcf < nbCF - 1; loopcf++) { rate0Ratio[loopcf] = (liborLMM[loopcf] + aLMM[loopcf]) / (liborLMM[loopcf] + 1 / deltaLMM[loopcf]); } for (int loopfact = 0; loopfact < nbFactor; loopfact++) { mu0[0][loopfact] = rate0Ratio[0] * gammaLMM[0][loopfact]; } for (int loopcf = 1; loopcf < nbCF - 1; loopcf++) { for (int loopfact = 0; loopfact < nbFactor; loopfact++) { mu0[loopcf][loopfact] = mu0[loopcf - 1][loopfact] + rate0Ratio[loopcf] * gammaLMM[loopcf][loopfact]; } } final double[] tau = new double[nbCF]; final double[] tau2 = new double[nbCF]; for (int loopcf = 0; loopcf < nbCF - 1; loopcf++) { for (int loopfact = 0; loopfact < nbFactor; loopfact++) { tau2[loopcf + 1] += mu0[loopcf][loopfact] * mu0[loopcf][loopfact]; } tau2[loopcf + 1] = tau2[loopcf + 1] * meanReversionImpact; tau[loopcf + 1] = Math.sqrt(tau2[loopcf + 1]); } double sumNum = -bM; double sumDen = 0; for (int loopcf = 0; loopcf < nbCF; loopcf++) { sumNum += dP[loopcf] - dP[loopcf] * tau2[loopcf] / 2.0; sumDen += dP[loopcf] * tau[loopcf]; } final double xBar = sumNum / sumDen; final double[] pM = new double[nbCF]; for (int loopcf = 0; loopcf < nbCF; loopcf++) { pM[loopcf] = p0[loopcf] * (1 - xBar * tau[loopcf] - tau2[loopcf] / 2.0); } final double[] liborM = new double[nbCF - 1]; final double[] alphaM = new double[nbCF]; for (int loopcf = 0; loopcf < nbCF - 1; loopcf++) { liborM[loopcf] = (pM[loopcf] / pM[loopcf + 1] - 1.0d) / deltaLMM[loopcf]; } for (int loopcf = 0; loopcf < nbCF; loopcf++) { alphaM[loopcf] = cfaMod[loopcf + 1] * pM[loopcf] / bM; } final double[] rateMRatio = new double[nbCF - 1]; final double[][] muM = new double[nbCF - 1][nbFactor]; for (int loopcf = 0; loopcf < nbCF - 1; loopcf++) { rateMRatio[loopcf] = (liborM[loopcf] + aLMM[loopcf]) / (liborM[loopcf] + 1 / deltaLMM[loopcf]); } for (int loopfact = 0; loopfact < nbFactor; loopfact++) { muM[0][loopfact] = rateMRatio[0] * gammaLMM[0][loopfact]; } for (int loopcf = 1; loopcf < nbCF - 1; loopcf++) { for (int loopfact = 0; loopfact < nbFactor; loopfact++) { muM[loopcf][loopfact] = muM[loopcf - 1][loopfact] + rateMRatio[loopcf] * gammaLMM[loopcf][loopfact]; } } double normSigmaM = 0; final double[] sigmaM = new double[nbFactor]; for (int loopfact = 0; loopfact < nbFactor; loopfact++) { for (int loopcf = 0; loopcf < nbCF - 1; loopcf++) { sigmaM[loopfact] += alphaM[loopcf + 1] * muM[loopcf][loopfact]; } normSigmaM += sigmaM[loopfact] * sigmaM[loopfact]; } final double impliedBlackVol = Math.sqrt(normSigmaM * meanReversionImpact); final EuropeanVanillaOption option = new EuropeanVanillaOption(bK, 1, isCall); final BlackPriceFunction blackFunction = new BlackPriceFunction(); final BlackFunctionData dataBlack = new BlackFunctionData(b0, 1.0, impliedBlackVol); final Function1D<BlackFunctionData, Double> func = blackFunction.getPriceFunction(option); final double pv = dfLMM[0] * func.evaluate(dataBlack); return MultipleCurrencyAmount.of(swaption.getUnderlyingSwap().getFirstLeg().getCurrency(), pv * (swaption.isLong() ? 1.0 : -1.0)); }
From source file:edu.cmu.tetrad.regression.LogisticRegression.java
/** * Regresses the single-column target onto the regressors which have been * previously set, generating a regression result. * <p>/*from www . j ava 2 s. c o m*/ * The target must be a two-valued variable with values 0 and 1. * <p> * This implements an iterative search. */ private Result regress(int[] target, String targetName, double[][] regressors, List<String> regressorNames) { double[][] x; double[] c1; int numRegressors = regressors.length; int numCases = target.length; // make a new matrix x with all the columns of regressors // but with first column all 1.0's. x = new double[numRegressors + 1][]; c1 = new double[numCases]; x[0] = c1; System.arraycopy(regressors, 0, x, 1, numRegressors); for (int i = 0; i < numCases; i++) { x[0][i] = 1.0; c1[i] = 1.0; } double[] xMeans = new double[numRegressors + 1]; double[] xStdDevs = new double[numRegressors + 1]; double[] y0 = new double[numCases]; double[] y1 = new double[numCases]; for (int i = 0; i < numCases; i++) { y0[i] = 0; y1[i] = 0; } int ny0 = 0; int ny1 = 0; int nc = 0; for (int i = 0; i < numCases; i++) { if (target[i] == 0.0) { y0[i] = 1; ny0++; } else { y1[i] = 1; ny1++; } nc += y0[i] + y1[i]; for (int j = 1; j <= numRegressors; j++) { xMeans[j] += (y0[i] + y1[i]) * x[j][i]; xStdDevs[j] += (y0[i] + y1[i]) * x[j][i] * x[j][i]; } } for (int j = 1; j <= numRegressors; j++) { xMeans[j] /= nc; xStdDevs[j] /= nc; xStdDevs[j] = Math.sqrt(Math.abs(xStdDevs[j] - xMeans[j] * xMeans[j])); } xMeans[0] = 0.0; xStdDevs[0] = 1.0; for (int i = 0; i < nc; i++) { for (int j = 1; j <= numRegressors; j++) { x[j][i] = (x[j][i] - xMeans[j]) / xStdDevs[j]; } } //report = report + ("Iteration history...\n"); double[] par = new double[numRegressors + 1]; double[] parStdErr = new double[numRegressors + 1]; double[] coefficients; par[0] = Math.log((double) ny1 / (double) ny0); for (int j = 1; j <= numRegressors; j++) { par[j] = 0.0; } double[][] arr = new double[numRegressors + 1][numRegressors + 2]; double lnV; double ln1mV; double llP = 2e+10; double ll = 1e+10; double llN = 0.0; while (Math.abs(llP - ll) > 1e-7) { /// 1e-7 llP = ll; ll = 0.0; for (int j = 0; j <= numRegressors; j++) { for (int k = j; k <= numRegressors + 1; k++) { arr[j][k] = 0.0; } } for (int i = 0; i < nc; i++) { double q; double v = par[0]; for (int j = 1; j <= numRegressors; j++) { v += par[j] * x[j][i]; } if (v > 15.0) { lnV = -Math.exp(-v); ln1mV = -v; q = Math.exp(-v); v = Math.exp(lnV); } else { if (v < -15.0) { lnV = v; ln1mV = -Math.exp(v); q = Math.exp(v); v = Math.exp(lnV); } else { v = 1.0 / (1 + Math.exp(-v)); lnV = Math.log(v); ln1mV = Math.log(1.0 - v); q = v * (1.0 - v); } } ll = ll - 2.0 * y1[i] * lnV - 2.0 * y0[i] * ln1mV; for (int j = 0; j <= numRegressors; j++) { double xij = x[j][i]; arr[j][numRegressors + 1] += xij * (y1[i] * (1.0 - v) + y0[i] * (-v)); for (int k = j; k <= numRegressors; k++) { arr[j][k] += xij * x[k][i] * q * (y0[i] + y1[i]); } } } if (llP == 1e+10) { llN = ll; } for (int j = 1; j <= numRegressors; j++) { for (int k = 0; k < j; k++) { arr[j][k] = arr[k][j]; } } for (int i = 0; i <= numRegressors; i++) { double s = arr[i][i]; arr[i][i] = 1.0; for (int k = 0; k <= numRegressors + 1; k++) { arr[i][k] = arr[i][k] / s; } for (int j = 0; j <= numRegressors; j++) { if (i != j) { s = arr[j][i]; arr[j][i] = 0.0; for (int k = 0; k <= numRegressors + 1; k++) { arr[j][k] = arr[j][k] - s * arr[i][k]; } } } } for (int j = 0; j <= numRegressors; j++) { par[j] += arr[j][numRegressors + 1]; } } //report = report + (" (Converged) \n"); // EdgeListGraph outgraph = new EdgeListGraph(); // Node targNode = new GraphNode(targetName); // outgraph.addNode(targNode); double chiSq = llN - ll; //Indicates whether each coefficient is significant at the alpha level. String[] sigMarker = new String[numRegressors]; double[] pValues = new double[numRegressors + 1]; double[] zScores = new double[numRegressors + 1]; for (int j = 1; j <= numRegressors; j++) { par[j] = par[j] / xStdDevs[j]; parStdErr[j] = Math.sqrt(arr[j][j]) / xStdDevs[j]; par[0] = par[0] - par[j] * xMeans[j]; double zScore = par[j] / parStdErr[j]; double prob = norm(Math.abs(zScore)); pValues[j] = prob; zScores[j] = zScore; } parStdErr[0] = Math.sqrt(arr[0][0]); double zScore = par[0] / parStdErr[0]; pValues[0] = norm(zScore); zScores[0] = zScore; double intercept = par[0]; coefficients = par; return new Result(targetName, regressorNames, xMeans, xStdDevs, numRegressors, ny0, ny1, coefficients, parStdErr, pValues, intercept, ll, sigMarker, chiSq, alpha); }
From source file:edu.duke.cs.osprey.voxq.QuadraticQFunction.java
private boolean normalizeDistr() { //want integral of exp(ax^2+bx+c) from xLo to xHi to equal 1 //adjust c accordingly. Return false if unsuccessful if (Math.abs(a) < 1e-14) {//linear case //for numerical reasons, we will treat this differently dependent on what b is if (Math.abs(b) < 1e-7) { //no noticeable variation across voxel: basically exp(c) * exp(b*(xHi+xLo)/2) * (xHi-xLo)=1 //(comes from linearization of cumulDistr; this form is actually correct to 2nd order in b c = -b * (xHi + xLo) / 2 - Math.log(xHi - xLo); } else if (b < 0) { //integral = exp(b*xLo+c) * (exp(b*(xHi-xLo))-1) / b double expl = b / (Math.exp(b * (xHi - xLo)) - 1); c = Math.log(expl) - b * xLo; } else {//b>=1e-7 //integral = exp(b*xHi+c) * (1-exp(b*(xLo-xHi))) / b double expl = b / (1 - Math.exp(b * (xLo - xHi))); c = Math.log(expl) - b * xHi; }/* ww w . j av a2 s . c o m*/ } else { double erfDiff = Erf.erf(cdErfArg(xLo), cdErfArg(xHi)); if (erfDiff != 0) {//erfDiff=0 prevents normalization; in this case will fail erfcInvNumericsTest double exponent = -Math.log(0.5 * erfDiff * Math.sqrt(-Math.PI / a)); c = exponent + 0.25 * b * b / a; } } double newNorm = cumulDistr(xHi); /* if( (Double.isInfinite(newNorm)||Double.isNaN(newNorm)) && Math.abs(a)>=1e-14 ) return false;//quadratic can't be normalized, so will be replaced by linear else if( Math.abs(newNorm-1) > 1e-5 ) throw new RuntimeException("ERROR: Unsuccessful normalization. a: "+a+" b: "+b+" c: "+c); else if( Double.isNaN(newNorm) ) throw new RuntimeException("ERROR: NaN normalization. a: "+a+" b: "+b+" c: "+c); */ if (Double.isNaN(newNorm) || Math.abs(newNorm - 1) > 1e-5) {//normalization failed if (Math.abs(a) >= 1e-14) {//this happens for nearly linear quadratic models //try explicitly linear instead return false; } else {//linear normalization failed...shouldn't happen for remotely realistic energies throw new RuntimeException("ERROR: Unsuccessful normalization. a: " + a + " b: " + b + " c: " + c + " newNorm (should be 1): " + newNorm); } } return true; }
From source file:de.tudarmstadt.ukp.dkpro.spelling.experiments.hoo2012.hoo2011.FixedCandidateTrigramProbabilityDetector.java
@Override protected double getSentenceProbability(List<String> words) throws AnalysisEngineProcessException { double sentenceProbability = 0.0; if (words.size() < 1) { return 0.0; }/*from w ww . j av a 2 s .c om*/ long nrOfUnigrams; try { nrOfUnigrams = provider.getNrOfTokens(); } catch (Exception e) { throw new AnalysisEngineProcessException(e); } List<String> trigrams = new ArrayList<String>(); // in the google n-grams this is not represented (only single BOS markers) // but I leave it in place in case we add another n-gram provider trigrams.add(NGramDetectorUtils.getTrigram(BOS, BOS, words.get(0))); if (words.size() > 1) { trigrams.add(NGramDetectorUtils.getTrigram(BOS, words.get(0), words.get(1))); } for (String trigram : new NGramStringIterable(words, 3, 3)) { trigrams.add(trigram); } // FIXME - implement backoff or linear interpolation for (String trigram : trigrams) { long trigramFreq = getNGramCount(trigram); String[] parts = StringUtils.split(trigram, " "); String bigram = StringUtils.join(Arrays.copyOfRange(parts, 0, 2), " "); long bigramFreq = getNGramCount(bigram); String unigram = StringUtils.join(Arrays.copyOfRange(parts, 0, 1), " "); long unigramFreq = getNGramCount(unigram); if (trigramFreq < 1) { trigramFreq = 1; } if (bigramFreq < 1) { bigramFreq = 1; } if (unigramFreq < 1) { unigramFreq = 1; } double trigramProb = Math.log((double) trigramFreq / bigramFreq); double bigramProb = Math.log((double) bigramFreq / unigramFreq); double unigramProb = Math.log((double) unigramFreq / nrOfUnigrams); double interpolated = (trigramProb + bigramProb + unigramProb) / 3.0; sentenceProbability += interpolated; } return Math.exp(sentenceProbability); }
From source file:edu.cmu.cs.lti.ark.fn.identification.latentmodel.TrainBatchModelDerThreaded.java
public void processBatch(int taskId, int start, int end) throws ExecutionException, IOException, ClassNotFoundException { int threadId = taskId % numThreads; logger.info("Processing batch:" + taskId + " thread ID:" + threadId); for (int targetIdx = start; targetIdx < end; targetIdx++) { int[][][] featureArray = getFeaturesForTarget(targetIdx); int featArrLen = featureArray.length; double exp[][] = new double[featArrLen][]; double sumExp[] = new double[featArrLen]; double totalExp = 0.0; for (int i = 0; i < featArrLen; i++) { exp[i] = new double[featureArray[i].length]; sumExp[i] = 0.0; // the partition function for (int j = 0; j < exp[i].length; j++) { double weiFeatSum = 0.0; int[] feats = featureArray[i][j]; for (int feat : feats) { weiFeatSum += params[feat]; }/*w w w.j ava2 s .c o m*/ exp[i][j] = Math.exp(weiFeatSum); sumExp[i] += exp[i][j]; } totalExp += sumExp[i]; } tValues[threadId] -= Math.log(sumExp[0]) - Math.log(totalExp); for (int i = 0; i < featArrLen; i++) { for (int j = 0; j < exp[i].length; j++) { double Y = 0.0; if (i == 0) { // the correct frame is always first Y = exp[i][j] / sumExp[i]; } double YMinusP = Y - (exp[i][j] / totalExp); int[] feats = featureArray[i][j]; for (int feat : feats) { tGradients[threadId][feat] -= YMinusP; } } } if (targetIdx % 100 == 0) { System.out.print("."); logger.info("" + targetIdx); } } if (useL2Regularization) { for (int i = 0; i < params.length; ++i) { final double weight = Math.log(params[i]); tValues[threadId] += lambda * (weight * weight); tGradients[threadId][i] += 2 * lambda * weight; } } }