List of usage examples for org.apache.commons.math3.linear ArrayRealVector ArrayRealVector
public ArrayRealVector(ArrayRealVector v) throws NullArgumentException
From source file:edu.utexas.cs.tactex.LWRCustomerPredictionTests.java
/** * testing LWR predictions - /*from ww w.j av a 2 s . co m*/ * based on matlab code for LWR * */ @Test public void testLWR() { // Ejml, with intercept SimpleMatrix X; SimpleMatrix y; SimpleMatrix x0; X = new SimpleMatrix(3, 2); y = new SimpleMatrix(3, 1); x0 = new SimpleMatrix(2, 1); // (with intercept) // initialize examples X.setColumn(0, 0, // X-column 1: (intercept) 1, 1, 1); X.setColumn(1, 0, // X-column 2: -1, 0, 1); y.setColumn(0, 0, // Y-column 3, 2, 4); x0.set(1); // initialize to 1s (will override non-intercept terms later) double tau = 0.5; // testing based on matlab results x0.set(1, 0, -1.5); double expected = 3.43855279053977; double actual = lwrCustNewEjml.LWRPredict(X, y, x0, tau); assertEquals("Matlab based expected value 1", expected, actual, 1e-6); x0.set(1, 0, -0.5); expected = 2.62107459906727; actual = lwrCustNewEjml.LWRPredict(X, y, x0, tau); assertEquals("Matlab based expected value 2", expected, actual, 1e-6); x0.set(1, 0, 0); expected = 2.63582467285126; actual = lwrCustNewEjml.LWRPredict(X, y, x0, tau); assertEquals("Matlab based expected value 3", expected, actual, 1e-6); x0.set(1, 0, 0.5); expected = 3.12107459906727; actual = lwrCustNewEjml.LWRPredict(X, y, x0, tau); assertEquals("Matlab based expected value 4", expected, actual, 1e-6); x0.set(1, 0, 1.5); expected = 4.93855279053977; actual = lwrCustNewEjml.LWRPredict(X, y, x0, tau); assertEquals("Matlab based expected value 5", expected, actual, 1e-6); // No intercept - should be same pred. between Ejml and Appache // build Ejml matrix with no intercept SimpleMatrix Xejml = X.extractMatrix(0, X.numRows(), 1, X.numCols()); // build appache matrix with no intercept ArrayRealVector Xappache = new ArrayRealVector(Xejml.numRows()); for (int i = 0; i < Xejml.numRows(); ++i) { Xappache.setEntry(i, Xejml.get(i, 0)); } ArrayRealVector yappache = new ArrayRealVector(y.numRows()); for (int i = 0; i < y.numRows(); ++i) { yappache.setEntry(i, y.get(i, 0)); } // initialize x0 x0 = new SimpleMatrix(1, 1); double x0val; double actualEjml; double actualAppache; x0val = -1.5; x0.set(0, 0, x0val); expected = 4.47403745685534; actualEjml = lwrCustNewEjml.LWRPredict(Xejml, y, x0, tau); actualAppache = lwrCustOldAppache.LWRPredict(Xappache, yappache, x0val, tau); assertEquals("NoIcpt - Matlab based - value 1, Ejml", expected, actualEjml, 1e-6); assertEquals("NoIcpt - Matlab based - value 1, Appache", expected, actualAppache, 1e-6); x0val = -0.5; x0.set(0, 0, x0val); expected = 1.08278977292259; actualEjml = lwrCustNewEjml.LWRPredict(Xejml, y, x0, tau); actualAppache = lwrCustOldAppache.LWRPredict(Xappache, yappache, x0val, tau); assertEquals("NoIcpt - Matlab based - value 2, Ejml", expected, actualEjml, 1e-6); assertEquals("NoIcpt - Matlab based - value 2, Appache", expected, actualAppache, 1e-6); x0val = 0; x0.set(0, 0, x0val); expected = 0; actualEjml = lwrCustNewEjml.LWRPredict(Xejml, y, x0, tau); actualAppache = lwrCustOldAppache.LWRPredict(Xappache, yappache, x0val, tau); assertEquals("NoIcpt - Matlab based - value 3, Ejml", expected, actualEjml, 1e-6); assertEquals("NoIcpt - Matlab based - value 3, Appache", expected, actualAppache, 1e-6); x0val = 0.5; x0.set(0, 0, x0val); expected = 1.58278977292259; actualEjml = lwrCustNewEjml.LWRPredict(Xejml, y, x0, tau); actualAppache = lwrCustOldAppache.LWRPredict(Xappache, yappache, x0val, tau); assertEquals("NoIcpt - Matlab based - value 4, Ejml", expected, actualEjml, 1e-6); assertEquals("NoIcpt - Matlab based - value 4, Appache", expected, actualAppache, 1e-6); x0val = 1.5; x0.set(0, 0, x0val); expected = 5.97403745685534; actualEjml = lwrCustNewEjml.LWRPredict(Xejml, y, x0, tau); actualAppache = lwrCustOldAppache.LWRPredict(Xappache, yappache, x0val, tau); assertEquals("NoIcpt - Matlab based - value 5, Ejml", expected, actualEjml, 1e-6); assertEquals("NoIcpt - Matlab based - value 5, Appache", expected, actualAppache, 1e-6); }
From source file:fp.overlapr.algorithmen.StressMajorization.java
/** * Fhrt die Stress-Majorization durch. siehe: Gansner, Koren, North: Graph * Drawing by Stress Majorization, 2004//from ww w .j a v a 2s.c o m * * @param graph * Graph, dessen Knoten-Positionen neu berechnet werden sollen * @param d * Matrix, die die idealen Distanzen (d_ij) zwischen den Knoten * enthlt * @return Matrix, die die neuen x- und y-Werte der einzelnen Knoten enthlt */ public static double[][] doStressMajorization(Graph graph, double[][] d) { int iter = 0; // X holen Array2DRowRealMatrix X = new Array2DRowRealMatrix(graph.getKnotenAnzahl(), 2); for (int i = 0; i < graph.getKnotenAnzahl(); i++) { X.setEntry(i, 0, graph.getKnoten().get(i).getX()); X.setEntry(i, 1, graph.getKnoten().get(i).getY()); } // D holen Array2DRowRealMatrix D = new Array2DRowRealMatrix(d); // W berechnen Array2DRowRealMatrix W = new Array2DRowRealMatrix(D.getRowDimension(), D.getColumnDimension()); W.walkInRowOrder(new DefaultRealMatrixChangingVisitor() { @Override public double visit(int row, int column, double value) { if (D.getEntry(row, column) == 0) { return 0.0; } else { return 1.0 / (D.getEntry(row, column) * D.getEntry(row, column)); } } }); // LW berechnen Array2DRowRealMatrix LW = new Array2DRowRealMatrix(D.getRowDimension(), D.getColumnDimension()); LW.walkInRowOrder(new DefaultRealMatrixChangingVisitor() { @Override public double visit(int row, int column, double value) { if (row != column) { return (-1) * W.getEntry(row, column); } else { return value; } } }); LW.walkInRowOrder(new DefaultRealMatrixChangingVisitor() { @Override public double visit(int row, int column, double value) { if (row == column) { double sum = 0; for (int k = 0; k < LW.getColumnDimension(); k++) { if (k != row) { sum = sum + W.getEntry(row, k); } } return sum; } else { return value; } } }); double[][] x = null; while (iter < ITER) { iter++; // LX berechnen Array2DRowRealMatrix LX = new Array2DRowRealMatrix(D.getRowDimension(), D.getColumnDimension()); LX.walkInRowOrder(new DefaultRealMatrixChangingVisitor() { @Override public double visit(int row, int column, double value) { if (row != column) { // norm 2 double term1 = FastMath.pow((X.getEntry(row, 0) - X.getEntry(column, 0)), 2); double term2 = FastMath.pow((X.getEntry(row, 1) - X.getEntry(column, 1)), 2); double abst = Math.sqrt(term1 + term2); return (-1) * W.getEntry(row, column) * D.getEntry(row, column) * inv(abst); } else { return value; } } }); LX.walkInRowOrder(new DefaultRealMatrixChangingVisitor() { @Override public double visit(int row, int column, double value) { if (row == column) { double sum = 0; for (int k = 0; k < LX.getColumnDimension(); k++) { if (k != row) { sum = sum + LX.getEntry(row, k); } } return (-1) * sum; } else { return value; } } }); /* * Lineare Gleichungssysteme lsen */ // x-Werte holen ArrayRealVector xWerte = new ArrayRealVector(X.getColumn(0)); // y-Werte holen ArrayRealVector yWerte = new ArrayRealVector(X.getColumn(1)); // b_x berechnen ArrayRealVector b_x = (ArrayRealVector) LX.operate(xWerte); // b_y berechnen ArrayRealVector b_y = (ArrayRealVector) LX.operate(yWerte); /* * CG-Verfahren anwenden */ // neue x-Werte berechnen mittels PCG // xWerte = conjugateGradientsMethod(LW, b_x, xWerte); // neue y-Werte berechnen mittels PCG // yWerte = conjugateGradientsMethod(LW, b_y, yWerte); ConjugateGradient cg = new ConjugateGradient(Integer.MAX_VALUE, TOL, false); xWerte = (ArrayRealVector) cg.solve(LW, JacobiPreconditioner.create(LW), b_x); yWerte = (ArrayRealVector) cg.solve(LW, JacobiPreconditioner.create(LW), b_y); /* * neue Positiones-Werte zurckgeben */ x = new double[X.getRowDimension()][2]; for (int i = 0; i < x.length; i++) { x[i][0] = xWerte.getEntry(i); x[i][1] = yWerte.getEntry(i); X.setEntry(i, 0, xWerte.getEntry(i)); X.setEntry(i, 1, yWerte.getEntry(i)); } } return x; }
From source file:info.rmarcus.birkhoffvonneumann.MatrixUtils.java
public static double[] normalize(double[] input) { ArrayRealVector arv = new ArrayRealVector(input); arv.mapDivideToSelf(arv.getNorm());// w w w.jav a 2 s .c o m return NullUtils.orThrow(arv.toArray(), () -> new BVNRuntimeException("Could not normalize array!")); }
From source file:gamlss.algorithm.Gamlss.java
/** This is to emulate the Gamlss algorithm where * the user specifies response variable vector * and design matrix.//w w w .ja v a 2s . co m * * @param y - vector of response variable values * @param designMatrices - design matrices for each * of the distribution parameters * @param smoothMatrices - smoother matrices for each * of the distribution parameters */ public Gamlss(final ArrayRealVector y, Hashtable<Integer, BlockRealMatrix> designMatrices, final HashMap<Integer, BlockRealMatrix> smoothMatrices) { GAMLSSFamilyDistribution distr = null; switch (DistributionSettings.DISTR) { case DistributionSettings.NO: distr = new NO(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.TF: distr = new TF(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.GA: distr = new GA(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.GT: distr = new GT(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.ST3: distr = new ST3(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.ST4: distr = new ST4(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.JSUo: distr = new JSUo(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.TF2: distr = new TF2(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.SST: distr = new SST(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.BCPE: distr = new BCPE(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.ST1: distr = new ST1(); distr.initialiseDistributionParameters(y); break; case DistributionSettings.PE: distr = new PE(); distr.initialiseDistributionParameters(y); break; default: System.err.println("The specific distribution " + "has not been implemented yet in Gamlss!"); } if (smoothMatrices != null) { Controls.SMOOTHING = true; } if (designMatrices == null) { designMatrices = new Hashtable<Integer, BlockRealMatrix>(); for (int i = 1; i < distr.getNumberOfDistribtionParameters() + 1; i++) { designMatrices.put(i, MatrixFunctions.buildInterceptMatrix(y.getDimension())); Controls.NO_INTERCEPT[i - 1] = true; } } ArrayRealVector w = new ArrayRealVector(y.getDimension()); for (int i = 0; i < w.getDimension(); i++) { w.setEntry(i, Controls.DEFAULT_WEIGHT); } switch (DistributionSettings.FITTING_ALG) { case DistributionSettings.RS: rs = new RSAlgorithm(distr, y, designMatrices, smoothMatrices, w); rs.functionRS(distr, y, w); break; case DistributionSettings.RS20CG20: rs = new RSAlgorithm(distr, y, designMatrices, smoothMatrices, w); rs.functionRS(distr, y, w); rs = null; cg = new CGAlgorithm(y.getDimension()); cg.setnCyc(Controls.GAMLSS_NUM_CYCLES); cg.CGfunction(distr, y, designMatrices, w); break; case DistributionSettings.CG: cg = new CGAlgorithm(y.getDimension()); cg.CGfunction(distr, y, designMatrices, w); break; case DistributionSettings.GO: break; default: System.err.println(" Cannot recognise the " + "fitting algorithm"); } int parNum = distr.getNumberOfDistribtionParameters(); resultsMatrix = new BlockRealMatrix(y.getDimension(), parNum); for (int i = 0; i < y.getDimension(); i++) { for (int j = 0; j < parNum; j++) { resultsMatrix.setEntry(i, j, distr.getDistributionParameter(j + 1).getEntry(i)); } } }
From source file:fi.smaa.jsmaa.model.PerCriterionMeasurementsTest.java
@Test public void testDeleteAlternatives() { ((MultivariateGaussianCriterionMeasurement) m.getCriterionMeasurement(criteria.get(0))) .setMeanVector(new ArrayRealVector(new double[] { 1.0, -1.0, -1.0 })); m.deleteAlternative(alternatives.get(0)); alternatives.remove(0);/*from www . j av a2 s.c o m*/ assertEquals(alternatives, m.getAlternatives()); assertEquals(alternatives, m.getCriterionMeasurement(criteria.get(0)).getAlternatives()); assertEquals(m.getCriterionMeasurement(criteria.get(0)).getRange(), ((ScaleCriterion) criteria.get(0)).getScale()); }
From source file:edu.stanford.cfuller.imageanalysistools.filter.MaximumSeparabilityThresholdingFilter.java
/** * Applies the filter to an Image, optionally turning on adaptive determination of the increment size used to find the threshold, and specifying a size for the threshold determination increment. * Turning on adaptive determination of the increment will generally make the threshold slightly less optimal, but can sometimes speed up the filtering, especially * for images with a large dynamic range. * <p>/* ww w . j a va 2 s . co m*/ * The increment size specified (in greylevels) will be used to determine the threshold only if adaptive determination is off; otherwise this parameter will be ignored. * @param im The Image to be thresholded; this will be overwritten by the thresholded Image. * @param adaptiveincrement true to turn on adaptive determination of the threshold increment; false to turn it off and use the default value * @param increment the increment size (in greylevels) to use for determining the threshold; must be positive. */ public void apply_ext(WritableImage im, boolean adaptiveincrement, int increment) { Histogram h = new Histogram(im); int thresholdValue = 0; final int numSteps = 1000; double best_eta = 0; int best_index = Integer.MAX_VALUE; int nonzerocounts = h.getTotalCounts() - h.getCounts(0); double meannonzero = h.getMeanNonzero(); ArrayRealVector omega_v = new ArrayRealVector(h.getMaxValue()); ArrayRealVector mu_v = new ArrayRealVector(h.getMaxValue()); int c = 0; if (adaptiveincrement) { increment = (int) ((h.getMaxValue() - h.getMinValue() + 1) * 1.0 / numSteps); if (increment < 1) increment = 1; } for (int k = h.getMinValue(); k < h.getMaxValue() + 1; k += increment) { if (k == 0) continue; omega_v.setEntry(c, h.getCumulativeCounts(k) * 1.0 / nonzerocounts); if (c == 0) { mu_v.setEntry(c, k * omega_v.getEntry(c)); } else { mu_v.setEntry(c, mu_v.getEntry(c - 1) + k * h.getCounts(k) * 1.0 / nonzerocounts); for (int i = k - increment + 1; i < k; i++) { mu_v.setEntry(c, mu_v.getEntry(c) + h.getCounts(i) * i * 1.0 / nonzerocounts); } } double omega = omega_v.getEntry(c); double mu = mu_v.getEntry(c); if (omega > 1e-8 && 1 - omega > 1e-8) { double eta = omega * (1 - omega) * Math.pow((meannonzero - mu) / (1 - omega) - mu / omega, 2); if (eta >= best_eta) { best_eta = eta; best_index = k; } //System.out.printf("%d, %f\n", k, eta); } c++; } thresholdValue = best_index; if (thresholdValue == Integer.MAX_VALUE) { thresholdValue = 0; } for (ImageCoordinate coord : im) { if (im.getValue(coord) < thresholdValue) im.setValue(coord, 0); } }
From source file:Gravitation_package.Gravitation02_ClassicWindowed.java
private void jButton1ActionPerformed(java.awt.event.ActionEvent evt) {//GEN-FIRST:event_jButton1ActionPerformed // TODO add your handling code here: SystemeDeParticules espace = new SystemeDeParticules(LARGEUR_ESPACE, HAUTEUR_ESPACE); espace.Creer(3);/*from ww w .ja va 2 s .c o m*/ final double x[] = { 0, 0 }; ArrayRealVector lieuTest = new ArrayRealVector(x); //Graphics g = new Graphics(): ArrayRealVector forceDeGravitation = espace.CalculerChampGravitationSurLieu(lieuTest); //System.out.println(forceDeGravitation.getX() + ";"+forceDeGravitation.getY()); }
From source file:iocomms.subpos.TrilaterationFunction.java
@Override public Pair<RealVector, RealMatrix> value(RealVector point) { // input/* ww w. j ava2 s .c om*/ double[] pointArray = point.toArray(); // output double[] resultPoint = new double[this.distances.length]; // compute least squares for (int i = 0; i < resultPoint.length; i++) { resultPoint[i] = 0.0; // calculate sum, add to overall for (int j = 0; j < pointArray.length; j++) { resultPoint[i] += (pointArray[j] - this.getPositions()[i][j]) * (pointArray[j] - this.getPositions()[i][j]); } resultPoint[i] -= (this.getDistances()[i]) * (this.getDistances()[i]); } RealMatrix jacobian = jacobian(point); return new Pair<RealVector, RealMatrix>(new ArrayRealVector(resultPoint), jacobian); }
From source file:edu.oregonstate.eecs.mcplan.ml.GaussianMixtureModel.java
private void init() { final int step = Math.max(1, n_ / k_); final double unif = 1.0 / k_; double acc = 0.0; final RandomPermutationIterator<double[]> r = new RandomPermutationIterator<double[]>(data_, rng_); final RandomPermutationIterator<double[]> rrepeat = new RandomPermutationIterator<double[]>(data_, r.permutation());//from w ww. j av a 2 s . c o m for (int i = 0; i < k_; ++i) { final RealVector mu = new ArrayRealVector(d_); for (int j = 0; j < step; ++j) { final double[] x = r.next(); final RealVector v = new ArrayRealVector(x); mu.combineToSelf(1.0, 1.0, v); } final double Zinv = 1.0 / step; mu.mapMultiplyToSelf(Zinv); RealMatrix Sigma = new Array2DRowRealMatrix(d_, d_); for (int j = 0; j < step; ++j) { final double[] x = rrepeat.next(); final RealVector v = new ArrayRealVector(x); v.combineToSelf(1.0, -1.0, mu); Sigma = Sigma.add(v.outerProduct(v)); } Sigma = Sigma.scalarMultiply(Zinv); pi_[i] = unif; acc += unif; mu_[i] = mu; Sigma_[i] = Sigma; //MatrixUtils.createRealIdentityMatrix( d_ ); } pi_[k_ - 1] += (1.0 - acc); // Round-off error }
From source file:com.joptimizer.algebra.CholeskySparseFactorizationTest.java
public void testSimple1_rescaling() throws Exception { log.debug("testSimple1_rescaling"); double[][] A = new double[][] { { 4, 0, 2, 2 }, { 0, 4, 2, -2 }, { 2, 2, 6, 0 }, { 2, -2, 0, 6 } }; //expected L/*from w w w. jav a 2s .co m*/ double[][] EL = new double[][] { { 2, 0, 0, 0 }, { 0, 2, 0, 0 }, { 1, 1, 2, 0 }, { 1, -1, 0, 2 } }; SparseDoubleMatrix2D ASparse = new SparseDoubleMatrix2D(A); CholeskySparseFactorization cs = new CholeskySparseFactorization(ASparse, new Matrix1NormRescaler()); cs.factorize(); DoubleMatrix2D L = cs.getL(); DoubleMatrix2D LT = cs.getLT(); log.debug("L : " + ArrayUtils.toString(L.toArray())); log.debug("LT: " + ArrayUtils.toString(LT.toArray())); RealMatrix ELMatrix = MatrixUtils.createRealMatrix(EL); RealMatrix LMatrix = MatrixUtils.createRealMatrix(L.toArray()); RealMatrix LTMatrix = MatrixUtils.createRealMatrix(LT.toArray()); assertEquals((ELMatrix.subtract(LMatrix).getNorm()), 0., Utils.getDoubleMachineEpsilon()); assertEquals((ELMatrix.subtract(LTMatrix.transpose()).getNorm()), 0., Utils.getDoubleMachineEpsilon()); //A.x = b double[] b = new double[] { 1, 2, 3, 4 }; double[] x = cs.solve(F1.make(b)).toArray(); //check the norm ||A.x-b|| double norm = new Array2DRowRealMatrix(A).operate(new ArrayRealVector(x)).subtract(new ArrayRealVector(b)) .getNorm(); log.debug("norm: " + norm); assertEquals(0., norm, 1.e-12); //check the scaled residual double residual = Utils.calculateScaledResidual(A, x, b); log.debug("residual: " + residual); assertEquals(0., residual, Utils.getDoubleMachineEpsilon()); }