List of usage examples for org.apache.commons.math3.linear ArrayRealVector getDimension
@Override public int getDimension()
From source file:gamlss.distributions.NO.java
/** Calculates initial value of sigma. * @param y - vector of values of response variable * @return vector of initial values of sigma */// w w w.j a v a 2s.c o m private ArrayRealVector setSigmaInitial(final ArrayRealVector y) { tempV = new ArrayRealVector(y.getDimension()); final double ySD = new StandardDeviation().evaluate(y.getDataRef()); tempV.set(ySD); return tempV; }
From source file:gamlss.utilities.MakeLinkFunction.java
/** * Calculates a linear predictor values of the distribution parameter vector according to the log link function * @param whichDistParameter - the distribution parameter * @return vector of linear predictor values *//*from ww w.ja v a2 s .co m*/ public ArrayRealVector log(ArrayRealVector whichDistParameter) { int size = whichDistParameter.getDimension(); double[] out = new double[size]; for (int i = 0; i < size; i++) { out[i] = FastMath.log(whichDistParameter.getEntry(i)); } return new ArrayRealVector(out, false); }
From source file:gamlss.utilities.MakeLinkFunction.java
/** * Calculates a linear predictor values of the distribution parameter vector according to the logshiftto2 link function * @param whichDistParameter - the distribution parameter * @return vector of linear predictor values *//*from w ww .ja va 2 s . c o m*/ public ArrayRealVector logShiftTo2(ArrayRealVector whichDistParameter) { int size = whichDistParameter.getDimension(); double[] out = new double[size]; for (int i = 0; i < size; i++) { //linkfun <- function(mu) log(mu-2+0.00001) # changed 6-7-12 out[i] = FastMath.log(whichDistParameter.getEntry(i) - 2 + 0.00001); } return new ArrayRealVector(out, false); }
From source file:nars.rl.util.ThreeDView.java
public ThreeDView(HyperassociativeMap h) { super();/*from w w w . j av a 2s .c o m*/ this.map = h; buildScene(); buildCamera(); buildAxes(); world.getChildren().addAll(space); double scale = 200.0; for (Object c : h.keys()) { ArrayRealVector v = h.getPosition((UIVertex) c); if (v.getDimension() >= 3) { double x = v.getEntry(0); double y = v.getEntry(1); double z = v.getEntry(2); addPoint(c, 40.0, x * scale, y * scale, z * scale); } } Scene scene = new Scene(root, 1024, 768, true); scene.setFill(Color.GREY); handleKeyboard(scene, world); handleMouse(scene, world); setScene(scene); scene.setCamera(camera); }
From source file:gamlss.distributions.GA.java
/** Calculates initial value of mu, by assumption these * values lie between observed data and the trend line. * @param y - vector of values of response variable * @return a vector of initial values of mu *///from w ww . ja v a 2 s. c o m private ArrayRealVector setMuInitial(final ArrayRealVector y) { //mu.initial = expression({mu <- (y+mean(y))/2}) size = y.getDimension(); double[] out = new double[size]; Mean mean = new Mean(); double yMean = mean.evaluate(y.getDataRef()); for (int i = 0; i < size; i++) { out[i] = (y.getEntry(i) + yMean) / 2; } return new ArrayRealVector(out, false); }
From source file:gamlss.utilities.MakeLinkFunction.java
/** * Calculates the values of distribution parameter Eta vector according to identity link function * @param eta - vector of linear predictor values * @return muEta vector/*from w ww .jav a 2 s.c om*/ */ public ArrayRealVector identityDistParameterEta(ArrayRealVector eta) { ArrayRealVector out = new ArrayRealVector(new double[eta.getDimension()], false); out.set(1); return out; }
From source file:gamlss.algorithm.GlimFitCG.java
public ArrayRealVector glimFitFunctionCG(GAMLSSFamilyDistribution distr, ArrayRealVector response, Hashtable<Integer, BlockRealMatrix> X, ArrayRealVector fv, ArrayRealVector weights, int whichDistParameter, double step, double offSet, double globalDeviance, Hashtable<String, ArrayRealVector> cgData, MakeLinkFunction makelink) { int itn = 0;//from w w w . j a v a 2 s . co m double[] temp = new double[0]; //G.dev.in <- G.dev+1 // double gDevIn = globalDeviance +1; //i.G.dev <- G.dev double gDev = globalDeviance; double gDevOld = globalDeviance + 1; //first.iter <- FALSE boolean firstIter = false; //while ( abs(G.dev.in -i.G.dev) > i.c.crit && i.iter < i.n.cyc ) while (FastMath.abs(gDevOld - gDev) > Controls.GLIM_CONV_CRIT && itn < Controls.GLIM_NUM_CYCLES) { //i.iter <- i.iter+1 itn = itn + 1; //for (int i = 1; i < distr.getNumberOfDistribtionParameters()+1; i++ ){ switch (whichDistParameter) { case DistributionSettings.MU: // whichDistParameter = DistributionSettings.MU; //adj.mu <- -(w.mu.sigma*(eta.sigma-eta.old.sigma)+w.mu.nu*(eta.nu-eta.old.nu)+w.mu.tau*(eta.tau-eta.old.tau))/w.mu adj = setAdj(cgData.get("wMuSigma"), cgData.get("eta" + DistributionSettings.SIGMA), cgData.get("etaOld" + DistributionSettings.SIGMA), cgData.get("wMuNu"), cgData.get("eta" + DistributionSettings.NU), cgData.get("etaOld" + DistributionSettings.NU), cgData.get("wMuTau"), cgData.get("eta" + DistributionSettings.TAU), cgData.get("etaOld" + DistributionSettings.TAU), cgData.get("w" + DistributionSettings.MU)); cgData.put("adj" + DistributionSettings.MU, adj); break; case DistributionSettings.SIGMA: // whichDistParameter = DistributionSettings.SIGMA; // adj.sigma <- -(w.mu.sigma*(eta.mu-eta.old.mu)+ w.sigma.nu*(eta.nu-eta.old.nu)+w.sigma.tau*(eta.tau-eta.old.tau))/w.sigma adj = setAdj(cgData.get("wMuSigma"), cgData.get("eta" + DistributionSettings.MU), cgData.get("etaOld" + DistributionSettings.MU), cgData.get("wSigmaNu"), cgData.get("eta" + DistributionSettings.NU), cgData.get("etaOld" + DistributionSettings.NU), cgData.get("wSigmaTau"), cgData.get("eta" + DistributionSettings.TAU), cgData.get("etaOld" + DistributionSettings.TAU), cgData.get("w" + DistributionSettings.SIGMA)); cgData.put("adj" + DistributionSettings.SIGMA, adj); break; case DistributionSettings.NU: // whichDistParameter = DistributionSettings.NU; break; case DistributionSettings.TAU: // whichDistParameter = DistributionSettings.TAU; break; } //wv.mu <- z.mu+adj.mu ArrayRealVector wv = setWV(cgData.get("z" + whichDistParameter), cgData.get("adj" + whichDistParameter)); //mu.fit <<- lm.wfit(x=mu.X,y=wv.mu,w=w.mu*w,method="qr") wls.newSampleData(wv, X.get(whichDistParameter).copy(), cgData.get("w" + whichDistParameter).ebeMultiply(weights)); ArrayRealVector fit = (ArrayRealVector) wls.calculateFittedValues(false); //System.out.println(wls.calculateBeta()); // bettas.put(whichDistParameter, (ArrayRealVector) wls.calculateBeta()); //mu.fit$eta <<- eta.mu <- mu.fit$fitted.values+mu.offset temp = new double[fit.getDimension()]; for (int j = 0; j < temp.length; j++) { temp[j] = fit.getEntry(j) + offSet; } //eta = new ArrayRealVector(temp,false); cgData.put("eta" + whichDistParameter, new ArrayRealVector(temp, false)); temp = null; //mu.fit$fv <<- mu <<- mu.object$linkinv(eta.mu) ArrayRealVector dParam = makelink.linkInv(distr.getDistributionParameterLink(whichDistParameter), cgData.get("eta" + whichDistParameter)); distr.setDistributionParameter(whichDistParameter, dParam); //mu.fit$wv <<- wv.mu //mu.fit$wt <<- w.mu //G.dev.in <- i.G.dev gDevOld = gDev; //G.dev.incr <- eval(G.dev.expr) //i.G.dev <- sum(w*G.dev.incr) gDev = weights.dotProduct(distr.globalDevianceIncreament(response)); if (gDev > gDevOld && itn >= 2 && Controls.AUTO_STEP) { for (int i = 0; i < 5; i++) { //eta.mu <- (eta.mu+eta.old.mu)/2 ArrayRealVector etaM = etaMean(cgData.get("eta" + whichDistParameter), cgData.get("etaOld" + whichDistParameter)); cgData.put("eta" + whichDistParameter, etaM); //mu <<- mu.object$linkinv(eta.mu) distr.setDistributionParameter(whichDistParameter, makelink.linkInv(distr.getDistributionParameterLink(whichDistParameter), etaM)); //if(length(who.mu) > 0) s.mu <- (s.mu+s.mu.old)/2 // } } } //if(i.trace) if (Controls.GLIM_TRACE) { //cat("CG inner iteration ", iter, ": Global Deviance = ",format(round(i.G.dev, 4)), " \n", sep = "") System.err.println("CG inner iteration " + itn + " : Global Deviance = " + gDev); } //if (i.G.dev > (G.dev.in+gd.tol) && iter >1 ) if (gDev > (gDevOld + Controls.GLOB_DEVIANCE_TOL) && itn > 1) { //stop(paste("The global deviance is increasing in the inner CG loop", "\n","Try different steps for the parameters or the model maybe inappropriate")) System.out.println( "The global deviance is increasing in the inner CG loop, try different steps for the parameters or the model maybe inappropriate"); break; } } //G.dev.old <- G.dev // gDevOld = gDev; //G.dev.incr <- eval(G.dev.expr) //G.dev <- sum(w*G.dev.incr) // gDev = weights.dotProduct(distr.globalDevianceIncreament(response)); //if (G.dev > G.dev.old && iter >= 2 && autostep == TRUE) return (ArrayRealVector) wls.calculateBeta(); }
From source file:edu.utexas.cs.tactex.utilityestimation.UtilityEstimatorDefaultForConsumption.java
/** * Core method for estimating utility//from ww w . j av a 2s. c o m */ public double estimateUtility(HashMap<TariffSpecification, HashMap<CustomerInfo, Integer>> tariffSubscriptions, HashMap<TariffSpecification, HashMap<CustomerInfo, Double>> predictedCustomerSubscriptions, HashMap<CustomerInfo, HashMap<TariffSpecification, Double>> customer2estimatedTariffCharges, HashMap<CustomerInfo, HashMap<TariffSpecification, ShiftedEnergyData>> customerTariff2ShiftedEnergy, HashMap<CustomerInfo, ArrayRealVector> customer2NonShiftedEnergy, MarketPredictionManager marketPredictionManager, CostCurvesPredictor costCurvesPredictor, int currentTimeslot/*, HashMap<CustomerInfo,ArrayRealVector> customer2estimatedEnergy*/) { // <= for print purposes log.debug( "estimateUtility(): currently assuming competing tariffs are not new and that my subscriptions are not going to change as a result of them"); // accumulate the final results //int expectedRecordLength = 7 * 24; // TODO use Timeservice here correctly double estTariffIncome = 0; double estTariffCosts = 0; double wholesaleCosts = 0; double balancingCosts = 0; double distributionCosts = 0; double withdrawCosts = 0; //double totalConsumption = 0; //double totalProduction = 0; final int predictionRecordLength = BrokerUtils.extractPredictionRecordLength(customerTariff2ShiftedEnergy); RealVector predictedEnergyRecord = new ArrayRealVector(predictionRecordLength); RealVector totalConsumptionEnergyRecord = new ArrayRealVector(predictionRecordLength); RealVector totalProductionEnergyRecord = new ArrayRealVector(predictionRecordLength); // for each customer get his usage prediction for (Entry<TariffSpecification, HashMap<CustomerInfo, Double>> entry : predictedCustomerSubscriptions .entrySet()) { TariffSpecification spec = entry.getKey(); HashMap<CustomerInfo, Double> subscriptions = entry.getValue(); for (Entry<CustomerInfo, Double> ce : subscriptions.entrySet()) { CustomerInfo customerInfo = ce.getKey(); Double subscribedPopulation = ce.getValue(); // Predicted total tariff cash flow. Sign is inverted since // evaluatedTariffs was computed from customers' perspective if (spec.getPowerType().isConsumption()) { estTariffIncome += -customer2estimatedTariffCharges.get(customerInfo).get(spec) * subscribedPopulation; } else if (spec.getPowerType().isProduction()) { estTariffCosts += -customer2estimatedTariffCharges.get(customerInfo).get(spec) * subscribedPopulation; } else { log.warn("Ignoring unknown powertype when computing tariffs income/costs: " + spec.getPowerType()); } // Predicted total energy RealVector energyPrediction = customerTariff2ShiftedEnergy.get(customerInfo).get(spec) .getShiftedEnergy().mapMultiply(subscribedPopulation); predictedEnergyRecord = predictedEnergyRecord.add(energyPrediction); // Predicted balancing cost balancingCosts += 0; log.debug("Ignoring balancing costs - assuming they are 0"); // Predicted withdraw costs (currently assuming everyone will pay). // sign is inverted since withdraw payment is from customer's perspective HashMap<CustomerInfo, Integer> cust2subs = tariffSubscriptions.get(spec); if (cust2subs != null) { Integer currentSubs = cust2subs.get(customerInfo); if (currentSubs != null && subscribedPopulation < currentSubs) { double withdraws = currentSubs - subscribedPopulation; withdrawCosts += -(withdraws * spec.getEarlyWithdrawPayment()); } } // Predicted total consumption and total production if (spec.getPowerType().isConsumption()) { totalConsumptionEnergyRecord = totalConsumptionEnergyRecord.add(energyPrediction); } else if (spec.getPowerType().isProduction()) { totalProductionEnergyRecord = totalProductionEnergyRecord.add(energyPrediction); } else { log.warn("Ignoring unknown powertype when computing distribution costs"); } } } // Predicted distribution costs log.debug("Ignoring balancing orders and curtailment when computing distribution costs"); distributionCosts = 0; double distributionFee = contextManager.getDistributionFee(); for (int i = 0; i < totalConsumptionEnergyRecord.getDimension(); ++i) { double totalTimeslotConsumption = Math.abs(totalConsumptionEnergyRecord.getEntry(i)); double totalTimeslotProduction = Math.abs(totalProductionEnergyRecord.getEntry(i)); distributionCosts += Math.max(totalTimeslotConsumption, totalTimeslotProduction) * distributionFee; } // Predicted wholesale costs (in one of the following two methods:) if (configuratorFactoryService.isUseCostCurves()) { // TODO: might not work for non-fixed rate competitor tariffs - // better to send a mapping of competitor tariffs inside // predictedCustomerSubscriptions, compute shifted predictions for them, // and use these predictions here. // compute energy of customers I don't have RealVector predictedCompetitorsEnergyRecord = new ArrayRealVector(predictionRecordLength); HashMap<CustomerInfo, HashMap<TariffSpecification, Double>> map = BrokerUtils .revertKeyMapping(predictedCustomerSubscriptions); for (CustomerInfo cust : map.keySet()) { double subsToOthers = cust.getPopulation() - BrokerUtils.sumMapValues(map.get(cust)); RealVector customerNonShiftedEnergy = customer2NonShiftedEnergy.get(cust).mapMultiply(subsToOthers); predictedCompetitorsEnergyRecord = predictedCompetitorsEnergyRecord.add(customerNonShiftedEnergy); } for (int i = 0; i < predictedEnergyRecord.getDimension(); ++i) { int futureTimeslot = currentTimeslot + i; double neededKwh = predictedEnergyRecord.getEntry(i); double competitorKwh = predictedCompetitorsEnergyRecord.getEntry(i); double unitCost = costCurvesPredictor.predictUnitCostKwh(currentTimeslot, futureTimeslot, neededKwh, competitorKwh); // NOTE: unitCost is signed (typically negative) wholesaleCosts += (unitCost + costCurvesPredictor.getFudgeFactorKwh(currentTimeslot)) * neededKwh; log.debug("cost-curve prediction: current " + currentTimeslot + " futureTimeslot " + futureTimeslot + " neededKwh " + neededKwh + " neededKwh + competitorKwh " + (neededKwh + competitorKwh) + " unitCost " + unitCost); } } else { // compute wholesale costs using this consumption and market manager marketprices prediction ArrayRealVector estimatedMarketPrices = marketPredictionManager.getPricePerKwhPredictionForAbout7Days(); // sanity check if (predictedEnergyRecord.getDimension() != estimatedMarketPrices.getDimension()) { log.error("Cannot compute utility - prediction periods of market and energy differ);"); return 0; } wholesaleCosts = -predictedEnergyRecord.dotProduct(estimatedMarketPrices); } log.info("estTariffIncome " + estTariffIncome); log.info("estTariffCosts " + estTariffCosts); log.info("wholesaleCosts " + wholesaleCosts); log.info("balancingCosts " + balancingCosts); log.info("distributionCosts " + distributionCosts); //log.info("old distributionCosts " + Math.max(totalProduction, totalConsumption) * contextManager.getDistributionFee()); log.info("withdrawCosts " + withdrawCosts); return estTariffIncome + estTariffCosts + wholesaleCosts + balancingCosts + distributionCosts + withdrawCosts; }
From source file:gamlss.distributions.GA.java
/** Calculates initial value of sigma. * @param y - vector of values of response variable * @return - a vector of initial values of sigma *//*from w w w . j a va 2 s . c o m*/ private ArrayRealVector setSigmaInitial(final ArrayRealVector y) { //sigma.initial = expression({sigma <- rep(1,length(y))}) , tempV = new ArrayRealVector(y.getDimension()); tempV.set(1); return tempV; }
From source file:gamlss.distributions.NO.java
/** Second cross derivative of likelihood function in * respect to mu and sigma (d2ldmdd = d2l/dmu*dsigma). * @param y - vector of values of response variable * @return a vector of Second cross derivative */// w w w . j a v a2 s . c o m private ArrayRealVector d2ldmdd(final ArrayRealVector y) { //d2ldmdd = function(y) rep(0,length(y) return new ArrayRealVector(y.getDimension()); // all elemnts are 0's); }