Estocásticos Programación con MCMC

Question

Mientras que hay muy superior a los métodos para la solución determinista LP problemas (por ejemplo, los algoritmos de punto interior), puede MCMC ser usado para resolver sus estocástico variantes?

Por estocásticos, me refiero a, por ejemplo, con las distribuciones de los coeficientes y / o limitaciones, y obtener una probabilidad posterior de la solución que representa:

• La incertidumbre en la solución dada la incertidumbre en los coeficientes/restricciones
• La incertidumbre en la solución teniendo en cuenta que MCMC le da una estimación

Por ejemplo:

\begin{equation*} \begin{array}{rrl} \mathbf{x}^* = \underset{\mathbf{x}}{\text{arg}\;\text{min}} & \mathbb{E}\,(c_1x_1 -3x_2)\\ \mbox{s.t.} & -x_1 +x_2 & \le b_1 \\ & x_1, x_2 & \geq 0 \end{array} \end{ecuación*}

y donde:

• $c_1 \sim N(2,0.5)$
• $b_1 \sim N(0,3)$

El objetivo es calcular los $pdf\mathbf(x^*) = \text{p}(x_1, x_2)$, reflejando nuestra incertidumbre sobre la optimalidad de la solución.

Estoy más familiarizado con PyMC, pero algo haría. Me inclino a creer que podíamos utilizar MRF potenciales para representar las restricciones y la función objetivo, pero no estoy seguro acerca de la manera de proceder para resolver el problema de optimización usando MCMC.

Answer 1

PyMC2 puede ser combinado con el LP solver de su elección para resolver estocástico LP problemas como este. Aquí está el código para hacerlo para este caso muy simple. He dejado una nota en cómo me iba a cambiar esto por una más compleja LP.

c1 = pm.Normal('c1', mu=2, tau=.5**-2)
c2 = -3
b1 = pm.Normal('b1', mu=0, tau=3.**-2)

@pm.deterministic
def x(c1=c1, c2=c2, b1=b1):
    # use an LP solver here for a complex problem
    arg_min = np.empty(2)
    min_val = np.inf
    for x1,x2 in [[0,0], [0, b1], [-b1, 0]]: # there are only three possible extreme points,
        if -x1 + x2 <= b1 and x1 >= 0 and x2 >= 0: # so check obj value at each valid one
            val = c1*x1 + c2*x2
            if val < min_val:
                min_val = val
                arg_min = [x1,x2]

    return np.array(arg_min, dtype=float)

Mira la extraña distribución conjunta de $(x_1, x_2)$:

Un cuaderno con todo el código para que esto está aquí.