close Warning: Can't synchronize with repository "(default)" (/var/svn/tolp does not appear to be a Subversion repository.). Look in the Trac log for more information.

Version 7 (modified by Víctor de Buen Remiro, 14 years ago) (diff)

--

  1. Generación de candidatos factibles en un politopo
    1. Descripción del problema
    2. Diseño del generador
      1. Definiciones
      2. Distribución de los candidatos
        1. Función de densidad
        2. Función generatriz

Generación de candidatos factibles en un politopo

Descripción del problema

Si tenemos que generar muestras que cumplan una serie de restricciones de desigualdad lineal

A x \ge a \wedge x\in\mathbb{R}^{n} \wedge a\in\mathbb{R}^{r} \wedge A\in\mathbb{R}^{r\times n}

caben dos posibilidades, utilizar un generador de candidatos que cumpla las restricciones por construcción o usar uno libre y luego rechazar los candidatos no factibles. La ventaja de este último es que puede ser simétrico y evita el cálculo de la verosimilitud del candidato pero el problema es que puede ser que tarde mucho en encontrar uno factible si el punto actual está demasiado cerca de la frontera, lo cual será muy habitual si el punto de máxima verosimilitud se encuentra fuera del politopo definido por las anteriores inecuaciones.

En la siguiente figura se observa un caso en el que el punto máximo verosímil sin restricciones (en verde) se encuentra fuera del politopo (en gris), por lo que el punto máximo-verosimil restringido (en rojo) se encuentra en la frontera del politopo. La cadena tenderá a estar cerca de ese punto y por tanto cerca de la frontera. Es evidente que cualquier entorno de generación simétrica de candidatos (en rosa), tendrá como máximo la mitad de área en la zona factible, lo cual es perfectamente asumible.

source:/tolp/OfficialTolArchiveNetwork/BysSampler/doc/image/RandWalk.InPolytope.chart.1026074815.png

Es claro que no es preciso que el punto se encuentre en la frontera para que sus entornos tengan incluso menos de la mitad, como ocurre con el entorno amarillo de la figura anterior. La pregunta es ¿qué pasaría si no estuviera en la frontera por inclumplir sólo una de las restricciones sino que incumpliera varias al mismo tiempo? Pues que la probabilidad de generar un punto factible decrecería exponencialmente, dependiendo también de los ángulos formados por los hiperplanos que definen cada restricción. Y eso ya no es asumible a nada que se tengan 3 ó mas restricciones incumplidas. En dos dimensiones sólo se pueden cruzar dos restricciones activas al mismo tiempo pero es fácil de extrapolar lo que ocurriría en espacios de dimensiones más altas:

source:/tolp/OfficialTolArchiveNetwork/BysSampler/doc/image/RandWalk.InPolytope.chart.726488582.png

Diseño del generador

Se hace necesario por lo tanto disponer de un generador de candidatos que por construcción estén siempre dentro del politopo, es decir un generador de candidatos factibles. Tal generador no podrá ser simétrico por lo que será necesario no sólo poder generar muestras sino también calcular su densidad de una forma eficiente.

Definiciones

Dado un punto x estrictamente interior al politopo

A x < a

la distancia al i-ésimo hiperplano viene dada por la fórmula

\left\langle x,\mathcal{H}\left(A_i,a_i\right)\right\rangle = \frac{\left|\left(\overset{n}{\underset{j=1}{\sum}}A_{ij}x_{j}\right)-a_{i}\right|}{\overset{n}{\underset{j=1}{\sum}}A_{ij}^{2}} \forall i=1\ldots r

Llamaremos distancia de un punto a la frontera del politopo a la menor de las distancias a cada uno de sus hiperplanos, la cual por ser el punto interior ha de ser forzosamente positiva

\left\langle x,\mathcal{P}\left(A,a\right)\right\rangle = \underset{i=1\ldots r}{min}\left\{ \left\langle x,\mathcal{H}\left(A_i,a_i\right)\right\rangle \} > 0

Distribución de los candidatos

Se propone utilizar un generador con distribución uniforme en una hiperesfera centrada en el último punto generado y que esté incluida estrictamente en el politopo lo cual será cierto si el radio es menor que la distancia del punto a la frontera:

\begin{array}{c} y=x + h u\\ A y \leq a\\ \left\Vert u\right\Vert = 1 \\ 0 < h \le \rho < \left\langle x,\mathcal{P}\left(A,a\right)\right\rangle \end{array}

Al igual que se hace con los generadores simétricos habituales es necesario fijar un tamaño de paso máximo s que para que se generen puntos no demasiado lejanos del acutual, de forma que el ratio de aceptación sea razonable. De esta forma se podría definir el radio de la hiperesfera como

\rho=\rho\left(x; s, \lambda\right)=min\left\{ s,\left\langle x, \lambda \mathcal{P}\left(A,a\right)\right\rangle \right\} \wedge 0<\lambda<1

Para poder muestrear puntos cercanos a la frontera sería conveniente que \lambda fuera cercano a 1, pero no tanto como para que los candidatos se quedaran demasiado cerca de esta oligado a que el radio se estrechara cada vez más.

Función de densidad

La densidad del generador será porporcional al inverso del volumen de la hiperesfera y su logaritmo será, salvo una constante

ln Q\left(x,y\right) = n ln \rho\left(x; s, \lambda\right)

el cual depende de la posición del punto de partida respecto a las fronteras del politopo y no es por tanto un generador simétrico necesariamente, auqnue sí puede serlo eventualmente si el la distancia del punto a la frontera es menor que el tamaño de paso actual.

Función generatriz

Para generar un candidato y a partir del actual x con esta distribución se procederá como sigue

  1. Primero se simulará un vector multinormal estandarizado

    w\sim Normal\left(0,I\right)\in\mathbb{R}^{n}

  2. Dividiendo ese vector por su norma euclídea se obtendrá una dirección unitaria uniforme en la frontera de la hiperesfera de radio 1 centrada en el origen

    u=\frac{w}{\left\Vert w\right\Vert _{2}}

  3. Luego se multiplica ese vector por un radio h con densidad directamente proporcional a la hipersuperficie de dimensión n-1 de la hiperesfera correspondiente

    \begin{array}{c} F\left(h\right)=\intop_{0}^{h}c\cdot t^{n-1}\mathbf{d}t=\left.\frac{c}{n}t^{n}\right]_{0}^{h}=\frac{c}{n}h^{n}\\ F\left(\rho\right)=\frac{c}{n}\rho^{n}=1\Rightarrow c=n\rho^{-n}\\ F\left(h\right)=\left(\frac{h}{\rho}\right)^{n}\end{array}