MOTA

Multi-objective tuning algorithm (MOTA)

class optTune.MOTA(objectiveFunctions, subproblems, CPV_lb, CPV_ub, CPV_validity_checks, sampleSizes, DE_F=2, DE_Cr=0.7, OFE_purtibation_factor=0.2, OFE_assessment_overshoot_function=<linearFunction y= 1.200000 * x + 100.000000 >, OFE_assessment_undershoot_function=<linearFunction y= 0.000000 * x + 0.000000 >, resampling_interruption_confidence=0.8, resampling_interruption_mode='reduce_max', boundConstrained=False, process_batch=<function _passfunction at 0x2b1ddd480230>, saveTo=None, saveInterval=-1, printFunction=<function to_stdout at 0x2b1de81526e0>, printLevel=2, record_X_hist=True, normalize_objective_values=True, post_iteration_function=<function _passfunction at 0x2b1ddd480230>, DE_F_vector_mutation=True, polynomial_similarity_mode=-1, simularity_exploitation_factor=2, simularity_fitting_threshold_ratio=0.2)

Multi-objective tuning algorithm

Required Args

• objectiveFunctions - contains the list of tuning objective functions. Each tuning objective (f) takes 3 arguments (CPV_array, assessment_OFE_budgets, randomSeed). f must returns two lists. The first list should return the utility measure such as the solution error (i.e f_found_opt - f_min) and should be decreasing, and the second list the number of objective function evaluations (OFEs) used in order to determine each element in the first list and should be increasing, and should if possible match the assessment_OFE_budgets (not required though). These lists can be thought of as the optimizer’s history. If an objective function has the addtoBatch attribute, each (CPV_array, assessment_OFE_budgets, randomSeed) about to be evaluated in passed to that function. Then the process batch_function is called, after which the objective function is called with the same input given to addtoBatch.
• subproblems - list of MOTA_sub_problems.
• CPV_lb - initialization lower bound for CPV tuning, i.e. numpy.array([CPV_1_init_lb, CPV_2_init_lb, ...])
• CPV_ub - numpy.array([CPV_1_init_ub, CPV_2_init_ub, ...])
• CPV_validity_checks - function used to check validity of candidate CPVs. Usage CPV_validity_checks(CPV_array, OFE_budget) returns tuple (valid, msg) where msg is string of why invalid. Should be a cheap function, as candidate CPVs are regenerated if they do not satisfy it. Use to prevent negative population sizes, populations size larger then OFE_budget checks, etcetera.
• sampleSizes - sample sizes used to generate and refined CPV utility values. For example if the sampleSizes are [5,15,30] all candidate CPVs will be sampled 5 times, the possible not dominated CPVs are then sampled another 15 times, and if still promising another 30 times. CPV which returned performances are therefore averaged over 50 independent instances.

Optional Args

• resampling_interruption_confidence - confidence level for interrupting the sample gathering process when performing re-sampling. Must be greater than 0.5 and less than or equal to 1
• resampling_interruption_mode - ‘reduce_max’ or ‘check_all’.
• OFE_purtubation_factor - when generating new candidate CPVs, OFEs close to the target OFE budget + 0.25 * rand_gausion* OFE_purtubation_factor * (log(max OFE budget)-log(min OFE budget)) are considered
• OFE_assessment_overshoot_function - when assessing a CPV tuple for a OFE budget of beta, this factor is used to control overshoot, beta_actual = OFE__assessment_undershot_function(beta)
• OFE_assessment_undershoot_function - similar to OFE__assessment_overshoot_function except controls minimum value
• DE_F - DE scaling factor
• DE_F_vector_mutation - change ‘DE_F*(x_1 - x_2)’ to ‘r()*DE_F*(x_1 - x_2)’ where r is vector consting randomly generated elements between 0 and 1, using a uniform distibution. Recommended, else diversity problem at start of MOTA optimization.
• DE_Cr - DE crossover probability
• boundConstrained - should CPV tuning search be bound constrainted between CPV_lb and CPV_ub
• process_batch - function is called after all the add_to_batch_functions have been called.
• saveTo - filename used to save optimization progress (one save per itteration), use None to disable saving
• saveInterval - number of iterations after which the progress of the tuning should be saved, negative numbers indicate do not saving during optimization.
• printLevel - 0 no printing, 1 only warnings, 2 overall info, 3 lots of info, 4 verbose (not implemented)
• printFunction - output is parsed into this function
• polynomial_similarity_mode - if an integer >= 1, then this CPV determines the order of the polynomial to be fitted for the similarity between polynomial fits of the paretoFronts (PFs) to construct the neighbourhood used for generating new candidate designs. A value of 0 indicates override_inter_PF_normalization, and use the update neighborhood for generating new candidate designs. A value of -1 indicates override_inter_PF_normalization, and every subproblem for generating new candidate designs. override_inter_PF_normalization - do not perform scaling correct based upon the fitted polynomials, during candidate CPV generation when using a CPV from a different PF to the target subproblem.
• simularity_exploitation_factor - controls how much information sharing take place between subproblems as a function of there similarity between their PFAs. Example values, 10 sharing only when very simular PFAs, 0 share equally regardless of simularity, -5 share with PFAs most dissimular. If function then, simularity_explotation_factor_iteration = simularity_exploitation_factor (subproblem.gamma / subproblem.gammaBudget)
• simularity_fitting_threshold_ratio - set the maximum scaling range, as (CPV_ub - CPV_lb)*simularity_scaling_threshold
• normalize_objective_values - normalize objective values so that utopia point is approximately 0 and nadir point is approximately 1 for all objectives
• post_iteration_function - at the end of each iteration this function is called with MOTA instance as the only arg.

MOTA subproblems

class optTune.MOTA_subproblem(w, target_OFE_budgets, gammaBudget, updateNeighbours, extra_termination_critea=[])

MOTA subproblem using Tchebycheff scalarization

Required Args * w - weights vector used for scalarization. * initial_target_OFE_budgets - OFE budgets at which tuning is desired to take place for the subproblem * gammaBudget - number of application layer evaluations to use on subproblem. i.e. if the algorithm being tuned is run 3 times with an OFE budget of 180 on 2 subproblem (w has 2 non-zero elements), then the gamma usage is 1080. * updateNeighbours - subproblems who PFAs are in the update neighbourhood when CPV tuples for this subproblem are evaluated

class optTune.MOTA_subproblem_weighted_sum(w, target_OFE_budgets, gammaBudget, updateNeighbours, extra_termination_critea=[])

MOTA subproblem using weighted sum scalarization

Required Args * w - weights vector used for scalarization. * initial_target_OFE_budgets - OFE budgets at which tuning is desired to take place for the subproblem * gammaBudget - number of application layer evaluations to use on subproblem. i.e. if the algorithm being tuned is run 3 times with an OFE budget of 180 on 2 subproblem (w has 2 non-zero elements), then the gamma usage is 1080. * updateNeighbours - subproblems who PFAs are in the update neighbourhood when CPV tuples for this subproblem are evaluated

optTune.generate_base_subproblem_list(n_f, target_OFE_budgets, gammaBudget, extra_termination_critea=, []subproblemClass=<class optTune.MOTA_code.subproblems.MOTA_subproblem at 0x2b1de8158258>)

Convenience function which return a subproblems for the weights [1,0, ..., 0], [0,1,...,0], ... [0,0,...,1] for the n_f objectives.

Example

Fortran code of optimization algorithm to tune under multiple OFE budgets

! written by Antoine Dymond, Jan 2012, based upon the Scipy.optimize.anneal code, the scipy license agreement is as follows:
!!$
!!$Copyright (c) 2001, 2002 Enthought, Inc.
!!$All rights reserved.
!!$
!!$Copyright (c) 2003-2009 SciPy Developers.
!!$All rights reserved.
!!$
!!$Redistribution and use in source and binary forms, with or without
!!$modification, are permitted provided that the following conditions are met:
!!$
!!$  a. Redistributions of source code must retain the above copyright notice,
!!$     this list of conditions and the following disclaimer.
!!$  b. Redistributions in binary form must reproduce the above copyright
!!$     notice, this list of conditions and the following disclaimer in the
!!$     documentation and/or other materials provided with the distribution.
!!$  c. Neither the name of the Enthought nor the names of its contributors
!!$     may be used to endorse or promote products derived from this software
!!$     without specific prior written permission.
!!$
!!$
!!$THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS "AS IS"
!!$AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
!!$IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE
!!$ARE DISCLAIMED. IN NO EVENT SHALL THE REGENTS OR CONTRIBUTORS BE LIABLE FOR
!!$ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL
!!$DAMAGES (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR
!!$SERVICES; LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
!!$CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
!!$LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY
!!$OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH
!!$DAMAGE.
!
! compile using
! $ f2py -c -m anneal_fortran anneal.f90


module anneal_module
  ! to compile using f2py, leave real type as default, needs to be compadiable between numpy and fortran.
  implicit None  ! fortran is case insensitive.
  !f2py forces all variable to lower case in linking module
  integer :: problem_id, objFun_evals
  real(8), allocatable, dimension(:) :: fval_hist, x_min
  integer, allocatable, dimension(:) :: eval_hist

contains

  function objFun(x, D)
    integer :: D, i
    real(8) :: x(D), objFun
    objFun = 0
    objFun_evals = objFun_evals + 1
    select case(problem_id)
    case (1) !general Rossenbrock
       do i = 1,D-1
          objFun = objFun + 100.0*(x(i+1) - x(i)**2.0)**2.0 + (1.0 - x(i))**2.0
       end do
    case (2) !sphere
       objFun = sum(x**2)
    end select
  end function objFun

  SUBROUTINE set_random_seed(seed_offset)
    !taken from the gfortran site. 
    !http://gcc.gnu.org/onlinedocs/gfortran/RANDOM_005fSEED.html#RANDOM_005fSEED
    INTEGER :: i, n, seed_offset
    INTEGER, DIMENSION(:), ALLOCATABLE :: seed   
    CALL RANDOM_SEED(size = n)
    ALLOCATE(seed(n))
    seed = seed_offset + 37 * (/ (i - 1, i = 1, n) /)
    CALL RANDOM_SEED(PUT = seed)
    DEALLOCATE(seed)
  END SUBROUTINE set_random_seed


  function rand_uniform()
    real(8) :: rand_uniform
    call random_number(rand_uniform)
  end function rand_uniform

  subroutine fast_sa_run(prob_id, x0, D, T0, dwell, m, n, quench, boltzmann, maxEvals, lower, upper, random_seed)
    integer :: prob_id, D, dwell, maxEvals, random_seed
    real(8) :: T0, m, n, quench, boltzmann
    real(8), dimension(D) :: x0, u, lower, upper 
    real(8) :: T, c, dE, y_j
    integer :: k, kMax, acceptTest, i, j
    real(8) :: current_state_fv, last_state_fv, best_state_fv
    real(8), dimension(D) ::current_state_xv, last_state_xv, best_state_xv


    call set_random_seed( random_seed )
    problem_id  = prob_id
    kMax = maxEvals/ dwell
    if (allocated(fval_hist)) deallocate(fval_hist, eval_hist )
    allocate(fval_hist(kMax), eval_hist(kMax))

    fval_hist = 0.0
    objFun_evals = 0
    T = T0
    c = m * exp(n * quench)
    last_state_xv = x0
    last_state_fv = objFun(x0, D)
    best_state_xv = x0
    best_state_fv = last_state_fv
    do k = 1,kMax
       do i = 1, dwell
          ! schedule.update_guess
          call random_number(u)
          do j = 1,D
             y_j = T * ( ( 1+1.0/T) ** abs(2*u(i)-1) - 1.0 )
             if ( (u(j) - 0.5) < 0 ) y_j = - y_j
             current_state_xv(j) = last_state_xv(j) + y_j*(upper(j) - lower(j))
             if (current_state_xv(j) < lower(j)) current_state_xv(j) = lower(j)
             if (current_state_xv(j) > upper(j)) current_state_xv(j) = upper(j)
          end do
          current_state_fv =  objFun(current_state_xv, D)
          dE = current_state_fv - last_state_fv
          !schedule.accept_test(dE)
          acceptTest = 0
          if (dE < 0) then
             acceptTest = 1
          else
             if  ( exp(-dE* 1.0 / boltzmann /T ) > rand_uniform() ) acceptTest = 1
          end if
          if ( acceptTest == 1) then
             last_state_xv = current_state_xv
             last_state_fv = current_state_fv
             if (last_state_fv < best_state_fv) then
                best_state_xv = last_state_xv
                best_state_fv = last_state_fv
             end if
          end if
       end do
       !tempreture update
       T = T0 * exp(-c*k*quench)
       !book keeping
       fval_hist(k) = best_state_fv
       eval_hist(k) = objFun_evals
    end do
    if (allocated(x_min)) deallocate(x_min)
    allocate(x_min(D))
    x_min = best_state_xv
  end subroutine fast_sa_run

end module anneal_module

Python code

import numpy, os
from optTune import MOTA,  get_F_vals_at_specified_OFE_budgets, linearFunction
from optTune import MOTA_subproblem_weighted_sum
print('Please note this example is writen for Linux, and requires gfortran')
if not os.path.exists('anneal_fortran.so'):
    os.system('f2py -c -m anneal_fortran anneal.f90')
from anneal_fortran import anneal_module
from matplotlib import pyplot

D = 12 #number of dimensions for Rosenbrock and Sphere problems
anneal_KWs = dict(
    t0 = 500.0,
    n = 1.0,
    quench = 1.0,
    boltzmann = 1.0,
)

def ros_ND( CPVs, OFE_budgets, randomSeed):
    anneal_module.fast_sa_run( 
        prob_id = 1,
        x0 = -2.048 + 2*2.048*numpy.random.rand(D),
        dwell = int(CPVs[0]), m = CPVs[1], 
        maxevals = max(OFE_budgets), 
        random_seed = randomSeed,
        lower = -2.048*numpy.ones(D),
        upper =  2.048*numpy.ones(D),
        **anneal_KWs )
    return get_F_vals_at_specified_OFE_budgets(F=anneal_module.fval_hist.copy(), E=anneal_module.eval_hist.copy(), E_desired=OFE_budgets)

def sphere_ND( CPVs, OFE_budgets, randomSeed):
    anneal_module.fast_sa_run( 
        prob_id = 2,
        x0 = -100 + 200*numpy.random.rand(D),
        dwell = int(CPVs[0]), m = CPVs[1], 
        maxevals = max(OFE_budgets), 
        random_seed = randomSeed,
        lower = -100.0*numpy.ones(D),
        upper =  100.0*numpy.ones(D),
        **anneal_KWs )
    return get_F_vals_at_specified_OFE_budgets(F=anneal_module.fval_hist.copy(), E=anneal_module.eval_hist.copy(), E_desired=OFE_budgets)

def CPV_valid(CPVs, OFE_budget):
    if CPVs[0] < 2:
        return False,'dwell,CPVs[0] < 2'
    if CPVs[1] < 0.0001:
        return False,'CPVs[1] < 0.0001'
    return True,''

subproblems =  []
for w_1 in numpy.linspace(0,1,5):
    subproblems.append( MOTA_subproblem_weighted_sum ( 
            w = numpy.array([w_1, 1-w_1]),
            target_OFE_budgets = numpy.logspace(1,3,30).astype(int),
            gammaBudget = 50*1000*50,
            updateNeighbours = [] #will do now
            ))
for i in range(len(subproblems)):
    subproblems[i].updateNeighbours = [subproblems[(i-1)%5],subproblems[(i+1)%5]]

tuningOpt = MOTA( 
    objectiveFunctions = [ros_ND, sphere_ND],
    subproblems = subproblems, 
    CPV_lb = numpy.array([30, 0.0]), 
    CPV_ub = numpy.array([50, 5.0]),
    CPV_validity_checks = CPV_valid,
    sampleSizes = [5,10,35], #resampling size of 30
    resampling_interruption_confidence = 0.6,
    printLevel=2,
    )
print(tuningOpt)

print('\nplotting optimal CPVs for different OFE budgets')

for j,S in enumerate(tuningOpt.subproblems):
    print(S)
    Fmin_values = [ d.fv[1] for d in S.PFA.designs ] 
    OFE_budgets = [ d.fv[0] for d in S.PFA.designs ]
    dwell_values =  [ int(d.xv[1]) for d in S.PFA.designs ]
    m_values = [ d.xv[2] for d in S.PFA.designs ]
    pyplot.subplot(1,5,j+1)

    p1 = pyplot.semilogx(OFE_budgets, dwell_values, 'g.')[0]
    pyplot.ylabel('dwell')
    pyplot.ylim( min(dwell_values) - 1, max( dwell_values) + 1)
    pyplot.twinx()
    p2 = pyplot.semilogx(OFE_budgets, m_values, 'bx')[0]
    pyplot.ylim( 0, max(m_values)*1.1)
    pyplot.ylabel('m (rate of cool)')
    pyplot.legend([p1,p2],['dwell','m'], loc='best')
    pyplot.xlabel('OFE budget')
    pyplot.xlim(min(OFE_budgets)-1,max(OFE_budgets)+60)
    pyplot.title('w = %s' % str(S.w))

pyplot.show()