/*
* This program is free software: you can redistribute it and/or modify
* it under the terms of the GNU General Public License as published by
* the Free Software Foundation, either version 3 of the License, or
* (at your option) any later version.
*
* This program is distributed in the hope that it will be useful,
* but WITHOUT ANY WARRANTY; without even the implied warranty of
* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
* GNU General Public License for more details.
*
* You should have received a copy of the GNU General Public License
* along with this program. If not, see <http://www.gnu.org/licenses/>.
*/
/*
* ConjugateGradientOptimization.java
* Copyright (C) 2012 University of Waikato, Hamilton, New Zealand
*
*/
package weka.core;
import weka.core.TechnicalInformation.Field;
import weka.core.TechnicalInformation.Type;
import java.util.Arrays;
/**
* This subclass of Optimization.java implements conjugate gradient
* descent rather than BFGS updates, by overriding findArgmin(), with
* the same tests for convergence, and applies the same line search
* code. Note that constraints are NOT actually supported. Using this
* class instead of Optimization.java can reduce runtime when there are
* many parameters.
*
* Uses the second hybrid method proposed in "An Efficient Hybrid
* Conjugate Gradient Method for Unconstrained Optimization" by Dai
* and Yuan (2001). See also information in the
* getTechnicalInformation() method.
*
* @author Eibe Frank
* @version $Revision: 8078 $
*/
public abstract class ConjugateGradientOptimization extends Optimization
implements RevisionHandler {
/**
* Returns an instance of a TechnicalInformation object, containing
* detailed information about the technical background of this class,
* e.g., paper reference or book this class is based on.
*
* @return the technical information about this class
*/
public TechnicalInformation getTechnicalInformation() {
TechnicalInformation result;
TechnicalInformation additional;
result = new TechnicalInformation(Type.ARTICLE);
result.setValue(Field.AUTHOR, "Y.H. Dai and Y. Yuan");
result.setValue(Field.YEAR, "2001");
result.setValue(Field.TITLE, "An Efficient Hybrid Conjugate Gradient Method for Unconstrained Optimization");
result.setValue(Field.JOURNAL, "Annals of Operations Research");
result.setValue(Field.VOLUME, "103");
result.setValue(Field.PAGES, "33-47");
additional = result.add(Type.ARTICLE);
result.setValue(Field.AUTHOR, "W.W. Hager and H. Zhang");
result.setValue(Field.YEAR, "2006");
result.setValue(Field.TITLE, "A survey of nonlinear conjugate gradient methods");
result.setValue(Field.JOURNAL, "Pacific Journal of Optimization");
result.setValue(Field.VOLUME, "2");
result.setValue(Field.PAGES, "35-58");
return result;
}
/**
* Constructor that sets MAXITS to 2000 by default and the parameter
* in the second weak Wolfe condition to 0.1.
*/
public ConjugateGradientOptimization() {
setMaxIteration(2000);
m_BETA = 0.1; // To make line search more exact, recommended for non-linear CGD
}
/**
* Main algorithm. NOTE: constraints are not actually supported.
*
* @param initX initial point of x, assuming no value's on the bound!
* @param constraints both arrays must contain Double.NaN
* @return the solution of x, null if number of iterations not enough
* @throws Exception if an error occurs
*/
public double[] findArgmin(double[] initX, double[][] constraints)
throws Exception{
int l = initX.length;
// Initial value of obj. function, gradient and inverse of the Hessian
m_f = objectiveFunction(initX);
if(Double.isNaN(m_f)) {
throw new Exception("Objective function value is NaN!");
}
// Get gradient at initial point
double[] grad = evaluateGradient(initX), oldGrad, oldX,
deltaGrad = new double[l], deltaX=new double[l],
direct = new double[l], x = new double[l];
// Turn gradient into direction and calculate squared length
double sum = 0;
for (int i = 0; i < grad.length; i++) {
direct[i] = -grad[i];
sum += grad[i] * grad[i];
}
// Same as in Optimization.java
double stpmax = m_STPMX * Math.max(Math.sqrt(sum), l);
boolean[] isFixed = new boolean[initX.length];
DynamicIntArray wsBdsIndx = new DynamicIntArray(initX.length);
double[][] consts = new double [2][initX.length];
for (int i = 0; i < initX.length; i++) {
if (!Double.isNaN(constraints[0][i]) || (!Double.isNaN(constraints[1][i]))) {
throw new Exception("Cannot deal with constraints, sorry.");
}
consts[0][i] = constraints[0][i];
consts[1][i] = constraints[1][i];
x[i] = initX[i];
}
boolean finished = false;
for (int step = 0; step < m_MAXITS; step++){
if (m_Debug) {
System.err.println("\nIteration # " + step + ":");
}
oldX = x;
oldGrad = grad;
// Make a copy of direction vector because it may get modified in lnsrch
double[] directB = Arrays.copyOf(direct, direct.length);
// Perform a line search based on new direction
m_IsZeroStep = false;
x = lnsrch(x, grad, directB, stpmax, isFixed, constraints, wsBdsIndx);
if (m_IsZeroStep) {
throw new Exception("Exiting due to zero step.");
}
// Check converge on x
boolean finish = false;
double test = 0.0;
for (int h = 0; h < x.length; h++){
deltaX[h] = x[h] - oldX[h];
double tmp = Math.abs(deltaX[h]) / Math.max(Math.abs(x[h]), 1.0);
if(tmp > test) test = tmp;
}
if (test < m_Zero){
if (m_Debug) {
System.err.println("\nDeltaX converged: " + test);
}
finished = true;
break;
}
// Check zero gradient
grad = evaluateGradient(x);
test = 0.0;
for (int g = 0; g < l; g++){
double tmp = Math.abs(grad[g]) * Math.max(Math.abs(directB[g]),1.0) / Math.max(Math.abs(m_f),1.0);
if (tmp > test) test = tmp;
}
if (test < m_Zero){
if (m_Debug) {
for (int i = 0; i < l; i++) {
System.out.println(grad[i] + " " + directB[i] + " " + m_f);
}
System.err.println("Gradient converged: " + test);
}
finished = true;
break;
}
// Calculate multiplier
double betaHSNumerator = 0, betaDYNumerator = 0;
double betaHSandDYDenominator = 0;
for (int i = 0; i < grad.length; i++) {
betaDYNumerator += grad[i] * grad[i];
betaHSNumerator += (grad[i] - oldGrad[i]) * grad[i];
betaHSandDYDenominator += (grad[i] - oldGrad[i]) * direct[i];
}
double betaHS = betaHSNumerator / betaHSandDYDenominator;
double betaDY = betaDYNumerator / betaHSandDYDenominator;
if (m_Debug) {
System.err.println("Beta HS: " + betaHS);
System.err.println("Beta DY: " + betaDY);
}
for (int i = 0; i < direct.length; i++) {
direct[i] = -grad[i] + Math.max(0, Math.min(betaHS, betaDY)) * direct[i];
}
}
if (finished) {
if (m_Debug){
System.err.println("Minimum found.");
}
m_f = objectiveFunction(x);
if(Double.isNaN(m_f)) {
throw new Exception("Objective function value is NaN!");
}
return x;
}
if(m_Debug) {
System.err.println("Cannot find minimum -- too many iterations!");
}
m_X = x;
return null;
}
}