jak zminimalizować funkcję z dyskretnymi wartościami zmiennych w scipy

Question

Próbuję zoptymalizować funkcję celu, która ma wiele zmiennych wejściowych (od 24 do 30). Zmienne te są próbkami trzech różnych zmiennych statystycznych, a wartości funkcji celu są wartościami prawdopodobieństwa testu t-Studenta. Funkcja błędu reprezentuje błąd (suma kwadratów różnic) między pożądanymi a rzeczywistymi prawdopodobieństwami testu t. Mogę zaakceptować rozwiązania, w których błąd jest mniejszy niż 1e-8, dla wszystkich trzech testów t.

Używałem scipy.optimize.fmin i działało świetnie. Istnieje wiele rozwiązań, w których funkcja docelowa stała się równa zero.

Problem polega na tym, że muszę znaleźć rozwiązanie, w którym zmienne mają wartość od 0 do 10,0 i są liczbami całkowitymi lub nie zawierają więcej niż jednej części ułamkowej. Przykłady poprawnych wartości to 0 10 3 5.5 6.8. Przykłady nieprawidłowych wartości: -3 2.23 30 lub 0.16666667.

Zdaję sobie sprawę, że istnieje przynajmniej jedno rozwiązanie, ponieważ wartości docelowe pochodzą z rzeczywistych zmierzonych danych. Pierwotne dane zostały utracone, a moim zadaniem jest je znaleźć. Ale nie wiem jak. Korzystanie z wersji próbnej/błędu nie jest opcją, ponieważ istnieje około 100 możliwych wartości dla każdej zmiennej, a biorąc pod uwagę liczbę zmiennych, liczba możliwych przypadków wynosi 100 ** 30, która jest zbyt duża. Używanie fmin jest świetne, ale nie działa z dyskretnymi wartościami.

Czy istnieje sposób rozwiązania tego problemu? Nie stanowi to problemu, jeśli potrzebuję uruchomić program przez wiele godzin, aby znaleźć rozwiązanie. Ale potrzebuję znaleźć rozwiązania dla około 10 wartości docelowych w ciągu kilku dni, a ja nie mam nowych pomysłów.

Oto przykład MWE:

import math 
import numpy 
import scipy.optimize 
import scipy.stats 
import sys 

def log(s): 
    sys.stdout.write(str(s)) 
    sys.stdout.flush() 

# List of target T values: TAB, TCA, TCB 
TARGETS = numpy.array([ 
    [0.05456834, 0.01510358, 0.15223353 ], # task 1 to solve 
    [0.15891875, 0.0083665,  0.00040262 ], # task 2 to solve 
]) 
MAX_ERR = 1e-10 # Maximum error in T values 
NMIN,NMAX = 8,10 # Number of samples for T probes. Inclusive. 

def fsq(x, t, n): 
    """Returns the differences between the target and the actual values.""" 
    a,b,c = x[0:n],x[n:2*n],x[2*n:3*n] 
    results = numpy.array([ 
     scipy.stats.ttest_rel(a,b)[1], # ab 
     scipy.stats.ttest_rel(c,a)[1], # ca 
     scipy.stats.ttest_rel(c,b)[1] # cb 
    ]) 
    # Sum of squares of diffs 
    return (results - t) 

def f(x, t, n): 
    """This is the target function that needs to be minimized.""" 
    return (fsq(x,t,n)**2).sum() 

def main(): 
    for tidx,t in enumerate(TARGETS): 
     print "=============================================" 
     print "Target %d/%d"%(tidx+1,len(TARGETS)) 
     for n in range(NMIN,NMAX+1): 
      log(" => n=%s "%n) 
      successful = False 
      tries = 0 
      factor = 0.1 
      while not successful: 
       x0 = numpy.random.random(3*n) * factor 
       x = scipy.optimize.fmin(f,x0, [t,n], xtol=MAX_ERR, ftol=MAX_ERR) 
       diffs = fsq(x,t,n) 
       successful = (numpy.abs(diffs)<MAX_ERR).all() 
       if successful: 
        log(" OK, error=[%s,%s,%s]\n"%(diffs[0],diffs[1],diffs[2])) 
        print " SOLUTION FOUND " 
        print x 
       else: 
        tries += 1 
        log(" FAILED, tries=%d\n"%tries) 
        print diffs 
        factor += 0.1 
        if tries>5: 
         print "!!!!!!!!!!!! GIVING UP !!!!!!!!!!!" 
         break 
if __name__ == "__main__": 
    main() 

Answer 1

Co próbujesz zrobić (jeśli zrozumiałem konfiguracji) nazywa programowanie całkowitą i jest NP-trudne; http://en.wikipedia.org/wiki/Integer_programming. Zdaję sobie sprawę, że nie szukasz rozwiązań całkowitych, ale jeśli pomnożysz wszystkie swoje dane wejściowe przez 10 i podzielisz docelową funkcję na 100, otrzymasz równorzędny problem, w którym dane wejściowe są liczbami całkowitymi. Chodzi o to, że twoje dane wejściowe są dyskretne.

Funkcja docelowa, z którą pracujesz, jest funkcją wypukłą, kwadratową i istnieją dobre algorytmy ograniczonej optymalizacji, które rozwiążą ją szybko dla wartościowań wejściowych w przedziale [0, 10]. Z tego możesz spróbować zaokrąglić lub sprawdzić wszystkie dopuszczalne punkty w pobliżu, ale jest ich 2^n, gdzie n to liczba wejść. Nawet jeśli to zrobisz, nie można zagwarantować, że optymalne rozwiązanie będzie jednym z tych punktów.

Istnieją algorytmy aproksymacyjne dla problemów z programowaniem całkowitoliczbowym i może się okazać, że czasami jedna z nich działa wystarczająco dobrze, aby dotrzeć do optymalnego punktu. Istnieje lista rzeczy, które możesz wypróbować w cytowanym przeze mnie artykule z Wikipedii, ale nie wiem, czy będziesz szczęśliwy, próbując rozwiązać ten problem.