Source code for schrodinger.math.multi_parameter_optimization

"""
Utility functions for MPO (Multi-Parameter Optimization). These provide a way to score
the desirability of some sort of object on a 0-1 scale based on individual
desirability scores of multiple properties of that object.

Use the `get_sigmoid` and `get_double_sigmoid` functions to get Sigmoid objects which
can be called with a parameter value to get the score. A list of scores and weights
can then be passed to get_weighted_score to compute an overall MPO score.
"""
import math
from collections import namedtuple

import numpy as np

# NOTE: The usage of the terms "good" and "bad" throughout are a useful way to
# discuss the different thresholds on a logistic curve in the context of MPO.
# But the fact that GOOD_Y > BAD_Y is ultimately arbitrary and these two
# values can equal anything between 0 and 1, exclusive. (GOOD_Y != BAD_Y)
GOOD_Y = .8  # the score of the threshold for good values
BAD_Y = .2  # the score of the threshold for bad values

# A rate and constant that determines a logistic sigmoid functions shape
Transition = namedtuple('Transition', ['center', 'rate'])
Transition.__new__.__defaults__ = (0.0, 1.0)


class Sigmoid:
    """
    A sigmoid transformation using a transition that determines the inflection
    point and the rate of change.

    Larger rate values lead to a faster transition.

    example usage:
    transition = Transition(2.5, 10.0)  # a steep 0 to 1 transition
    transform = Sigmoid(transition)
    transformed_valued = transform(original_value)

    :ivar transition: the transition
    :vartype transition: Transition

    ::

        rate positive                rate negative
        1      ---                   1 ---
              /                           \
        0  ---                       0     ---
           0 1 2 3 4...                0 1 2 3 4...
    """

    def __init__(self, transition):
        self._transition = transition

    def __call__(self, value):
        center, rate = self._transition
        return 1 / (1 + np.exp(-rate * (value - center)))


class DoubleSigmoid:
    """
    A double sigmoid transformation using a left and right transition to
    determine the inflection points and the rate of change. The left and
    right transitions are expected to have opposite signs for the rates.

    Larger rate values lead to a faster transition.

    example usage:
    a_side = Transform(1.5, 50.0)  # a very steep 0 to 1 transition
    b_side = Transition(10.0, -1.0) # a regular 1 to 0  transition
    mpo = DoubleSigmoid(a_side, b_side)
    transformed_valued = mpo(original_value)

    :ivar transition_a: the 'left-hand side' transition
    :vartype transition_a: Transition
    :ivar transition_b: the 'right-hand side' transition
    :vartype transition_b: Transition
    """

    def __init__(self, transition_a, transition_b):
        if np.sign(transition_a.rate) == np.sign(transition_b.rate):
            raise ValueError("Sigmoid curve rates must have opposite sign")
        self._sigmoid_a = Sigmoid(transition_a)
        self._sigmoid_b = Sigmoid(transition_b)
        self._intersection = get_intersection(transition_a, transition_b)

    def __call__(self, value):
        sigmoid = self._sigmoid_a if value < self._intersection else self._sigmoid_b
        return sigmoid(value)


def get_intersection(transition_a, transition_b):
    """
    Get the intersection point of two logistic functions given their center
    (inflection) points and rate constants.

    !!! Rate of transition_a and transition_b cannot be equal or a
    ZeroDivisionError will be raised. (Two curves with same rates will never
    intersect unless their centers are equal)

    :type transition_a: Transition
    :type transition_b: Transition
    :rtype: float
    """
    c0, r0 = transition_a
    c1, r1 = transition_b
    return (r1 * c1 - r0 * c0) / (r1 - r0)


def get_rate(good, bad):
    """
    Get the rate for the logistic function given the x values corresponding to the
    "good" and "bad" thresholds.

    :param good: The x value corresponding to the good threshold (`GOOD_Y`)
    :type good: float
    :param bad: The x value corresponding to the bad threshold (`BAD_Y`)
    :type bad: float
    :return: The rate
    :rtype: float
    """
    return (-math.log((1 / BAD_Y - 1) / (1 / GOOD_Y - 1))) / (bad - good)


def get_center(good, bad):
    """
    Get the center (inflection point) for the logistic function given the x
    values corresponding to the "good" and "bad" thresholds.

    :param good: The x value corresponding to the good threshold (`GOOD_Y`)
    :type good: float
    :param bad: The x value corresponding to the bad threshold (`BAD_Y`)
    :type bad: float
    :return: The center point
    :rtype: float
    """
    l_bad = math.log(1 / BAD_Y - 1)
    l_good = math.log(1 / GOOD_Y - 1)
    return (l_bad * good - l_good * bad) / (l_bad - l_good)


def get_sigmoid(good, bad):
    """
    Get a sigmoid logistic function using the given "good" and "bad" cutoffs

    :param good: The x value corresponding to the good threshold (`GOOD_Y`)
    :type good: float
    :param bad: The x value corresponding to the bad threshold (`BAD_Y`)
    :type bad: float
    :return: a sigmoid function which can be called with an x value to get
        the desirability
    :rtype: Sigmoid
    """
    return Sigmoid(Transition(get_center(good, bad), get_rate(good, bad)))


def get_double_sigmoid(good1, bad1, good2, bad2):
    """
    Get a double sigmoid function comprised of two logistic functions that
    fit two points on either function.

    :param good1: The x value corresponding to the good threshold (`GOOD_Y`) of the first sigmoid
    :type good1: float
    :param bad1: The x value corresponding to the bad threshold (`BAD_Y`) of the first sigmoid
    :type bad1: float
    :param good2: The x value corresponding to the good threshold (`GOOD_Y`) of the second sigmoid
    :type good2: float
    :param bad2: The x value corresponding to the bad threshold (`BAD_Y`) of the second sigmoid
    :type bad2: float
    :return: the double sigmoid function which can be called with an x value
        to get the desirability
    :rtype: DoubleSigmoid
    """
    if (bad1 < good1) is (bad2 < good2):
        raise ValueError("Sigmoid curves must have opposite sign")
    transition_a = Transition(get_center(good1, bad1), get_rate(good1, bad1))
    transition_b = Transition(get_center(good2, bad2), get_rate(good2, bad2))
    return DoubleSigmoid(transition_a, transition_b)


def get_weighted_score(scores, weights):
    """
    Return the weighted geometric mean of the given score

    :type scores: list[float]
    :type weights: list[float]
    :rtype: float
    """
    total_weight = sum(weights)
    if total_weight == 0:
        return 0.
    product = 1
    for score, weight in zip(scores, weights):
        product *= pow(score, weight)
    return pow(product, 1 / total_weight)