не сходимость мистического оптимизатора - PullRequest
0 голосов
/ 04 мая 2019

Я пытаюсь оптимизировать функцию objective.

full_constr_data состоит из 6 типов целей, каждая цель делится на годы, каждый год представлен данными на основе проекта.Итак, я взвешиваю full_constr_data по проекту с помощью аргумента x функции objective, например, full_constr_data[0][2][3] * x[3] означает, что "цель # 0, год # 2, проект # 3 взвешена на x[3] Результатыхранится в переменной full_constr_data_weighted.

Следующим шагом является суммирование по проектам full_constr_data_weighted. Пример, суммирование всех проектов для цели № 0, год № 2:

full_constr_data_weighted[0][2][0] + full_constr_data_weighted[0][2][1] + ... + full_constr_data_weighted[0][2][n]

где 'n' - общее количество проектов. Данные хранятся в переменной full_sum.

После этого я вычисляю вероятности. Я беру квантили из переменной constr_mod и основываюПо значению рассчитайте вероятности превышения квантиля. constr_mod и full_sum имеют абсолютно одинаковую структуру, однако для каждой цели # и следующего года # constr_mode содержит одно значение, тогда как full_sum имеет вектор значений (распределение). Вычисленные вероятности сохраняются в переменной my_prob.

Наконец, я суммирую все вероятности в my_prob. Эту сумму я должен оптимизировать: сделать ее как можно большей (обратите внимание назнак минус в выражении return).

ОптимВ задаче iizing есть одно ограничение неравенства: сумма вектора obj, взвешенная на x, должна быть больше 1000. Интерпретировать obj как NPV для каждого проекта.

Я использую diffev2 из пакетаmystic.

Переменные full_constr_data, constr_mod, pen_mult, obj хранятся в my_data.spydata file: загрузка через Google Drive

К сожалению, оптимизация не сходится:

Generation 11980 has ChiSquare: inf
Generation 11990 has ChiSquare: inf
Generation 12000 has ChiSquare: inf
STOP("EvaluationLimits with {'evaluations': 1200000, 'generations': 12000}")

Есть предложения, как решить эту невыпуклую проблему?

import numpy as np
from mystic.solvers import diffev2
from mystic.monitors import VerboseMonitor
import mystic.symbolic as ms

def objective(x):
    # 'full_constr_data' weighted by argument 'x'
    full_constr_data_weighted = []
    for i in range(len(full_constr_data)):
        temp = []
        for k in range(len(full_constr_data[i])):
            temp.append( [ full_constr_data[i][k][m] * x[m] \
                          for m in range(len(full_constr_data[i][k])) ] )
        full_constr_data_weighted.append(temp)

    # Project-wise sum of weighted data
    full_sum = []
    for i in range(len(full_constr_data_weighted)):
        temp = []
        for j in range(len(full_constr_data_weighted[i])):
            temp2 = np.array( [ 0 for m in range(len(full_constr_data_weighted[i][j][0])) ] )
            for k in range(len(full_constr_data_weighted[i][j])):
                temp2 = temp2 + full_constr_data_weighted[i][j][k]
            temp.append(temp2)
        full_sum.append(temp)

    # Probability calculation
    my_prob = []
    for i in range(len(full_sum)):
        temp = []
        for j in range(len(full_sum[i])):
            temp.append(sum(full_sum[i][j] > constr_mod[i][j]) / len(full_sum[i][j]))
        my_prob.append(np.array(temp))

    # Probability data weighted by 'pen_mult'
    my_prob_uweighted = list(np.array(my_prob) * np.array(pen_mult))

    # Sum of all weighted probability data (function to maximize)
    sum_prob = sum([sum(my_prob_uweighted[i]) for i in range(len(my_prob_uweighted))])

    return -sum_prob



# Inequality constraint
equation = 'x0*{0} + x1*{1} + x2*{2} + x3*{3} + x4*{4} + x5*{5} + x6*{6} + x7*{7} + x8*{8} + x9*{9} + x10*{10} + x11*{11} + x12*{12} + x13*{13} + x14*{14} + x15*{15} + x16*{16} + x17*{17} + x18*{18} + x19*{19} + x20*{20} + x21*{21} + x22*{22} + x23*{23} + x24*{24} + x25*{25} + x26*{26} + x27*{27} + x28*{28} + x29*{29} >= {30}'\
.format(obj[0], obj[1], obj[2], obj[3], obj[4], obj[5], obj[6], obj[7], obj[8], obj[9],
        obj[10], obj[11], obj[12], obj[13], obj[14], obj[15], obj[16], obj[17], obj[18], obj[19],
        obj[20], obj[21], obj[22], obj[23], obj[24], obj[25], obj[26], obj[27], obj[28], obj[29],
        1000)


cf = ms.generate_constraint(ms.generate_solvers(ms.simplify(equation)))

bounds = [(0,1)]*30

mon = VerboseMonitor(10)

result = diffev2(objective, x0=bounds, bounds=bounds, constraints=cf, \
                 npop=40, gtol=200, disp=False, full_output=True, itermon=mon)
...