Aproximación estocástica

Familia de métodos iterativos

Los métodos de aproximación estocástica son una familia de métodos iterativos que se utilizan normalmente para problemas de búsqueda de raíces o de optimización . Las reglas de actualización recursiva de los métodos de aproximación estocástica se pueden utilizar, entre otras cosas, para resolver sistemas lineales cuando los datos recopilados están alterados por el ruido, o para aproximar valores extremos de funciones que no se pueden calcular directamente, sino que solo se pueden estimar mediante observaciones ruidosas.

En pocas palabras, los algoritmos de aproximación estocástica tratan con una función de la forma que es el valor esperado de una función que depende de una variable aleatoria . El objetivo es recuperar propiedades de dicha función sin evaluarla directamente. En cambio, los algoritmos de aproximación estocástica utilizan muestras aleatorias de para aproximar de manera eficiente propiedades de como ceros o extremos. ${\textstyle f(\theta )=\operatorname {E} _{\xi }[F(\theta ,\xi )]}$ ${\textstyle \xi }$ ${\textstyle f}$ ${\textstyle F(\theta ,\xi )}$ ${\textstyle f}$

Recientemente, las aproximaciones estocásticas han encontrado amplias aplicaciones en los campos de la estadística y el aprendizaje automático, especialmente en entornos con big data . Estas aplicaciones van desde métodos y algoritmos de optimización estocástica , hasta formas en línea del algoritmo EM , aprendizaje de refuerzo mediante diferencias temporales y aprendizaje profundo , entre otros. ^[1] Los algoritmos de aproximación estocástica también se han utilizado en las ciencias sociales para describir dinámicas colectivas: el juego ficticio en la teoría del aprendizaje y los algoritmos de consenso se pueden estudiar utilizando su teoría. ^[2]

Los primeros algoritmos prototípicos de este tipo son los algoritmos Robbins-Monro y Kiefer-Wolfowitz , introducidos respectivamente en 1951 y 1952.

Algoritmo de Robbins-Monro

El algoritmo Robbins–Monro, introducido en 1951 por Herbert Robbins y Sutton Monro ^[3] , presentó una metodología para resolver un problema de búsqueda de raíces, donde la función se representa como un valor esperado. Supongamos que tenemos una función y una constante , de modo que la ecuación tiene una raíz única en . Se supone que, si bien no podemos observar directamente la función , podemos obtener mediciones de la variable aleatoria donde . La estructura del algoritmo es generar iteraciones de la forma: ${\textstyle M(\theta )}$ ${\textstyle \alpha }$ ${\textstyle M(\theta )=\alpha }$ ${\textstyle \theta ^{*}}$ ${\textstyle M(\theta )}$ ${\textstyle N(\theta )}$ ${\textstyle \operatorname {E} [N(\theta )]=M(\theta )}$

${\displaystyle \theta _{n+1}=\theta _{n}-a_{n}(N(\theta _{n})-\alpha )}$

Aquí, hay una secuencia de tamaños de paso positivos. Robbins y Monro demostraron ^{[3] , el Teorema 2} que converge en (y por lo tanto también en probabilidad) a , y Blum ^[4] demostró más tarde que la convergencia es en realidad con probabilidad uno, siempre que: ${\displaystyle a_{1},a_{2},\dots }$ ${\displaystyle \theta _{n}}$ ${\displaystyle L^{2}}$ ${\displaystyle \theta ^{*}}$

• ${\textstyle N(\theta )}$está uniformemente acotado,
• ${\textstyle M(\theta )}$no es decreciente,
• ${\textstyle M'(\theta ^{*})}$existe y es positivo, y
• La secuencia satisface los siguientes requisitos: ${\textstyle a_{n}}$

${\displaystyle \qquad \sum _{n=0}^{\infty }a_{n}=\infty \quad {\mbox{ and }}\quad \sum _{n=0}^{\infty }a_{n}^{2}<\infty \quad }$

Una secuencia particular de pasos que satisface estas condiciones, y fue sugerida por Robbins–Monro, tiene la forma: , para . Otras series son posibles pero para promediar el ruido en , se debe cumplir la condición anterior. ${\textstyle a_{n}=a/n}$ ${\textstyle a>0}$ ${\textstyle N(\theta )}$

Resultados de complejidad

1. Si es dos veces continuamente diferenciable, y fuertemente convexo, y el minimizador de pertenece al interior de , entonces el algoritmo de Robbins-Monro alcanzará la tasa de convergencia asintóticamente óptima, con respecto a la función objetivo, siendo , donde es el valor mínimo de sobre . ^{[5] [6]} ${\textstyle f(\theta )}$ ${\textstyle f(\theta )}$ ${\textstyle \Theta }$ ${\textstyle \operatorname {E} [f(\theta _{n})-f^{*}]=O(1/n)}$ ${\textstyle f^{*}}$ ${\textstyle f(\theta )}$ ${\textstyle \theta \in \Theta }$
2. Por el contrario, en el caso convexo general, donde carecemos tanto del supuesto de suavidad como de convexidad fuerte, Nemirovski y Yudin ^[7] han demostrado que la tasa de convergencia asintóticamente óptima, con respecto a los valores de la función objetivo, es . También han demostrado que esta tasa no se puede mejorar. ${\textstyle O(1/{\sqrt {n}})}$

Desarrollos posteriores y promedio de Polyak-Ruppert

Si bien el algoritmo Robbins-Monro es teóricamente capaz de lograr resultados bajo el supuesto de una diferenciabilidad continua doble y una fuerte convexidad, puede tener un desempeño bastante deficiente en la implementación. Esto se debe principalmente al hecho de que el algoritmo es muy sensible a la elección de la secuencia de tamaño de paso, y la supuesta política de tamaño de paso asintóticamente óptima puede ser bastante perjudicial al principio. ^{[6] [8]} ${\textstyle O(1/n)}$

Chung (1954) ^[9] y Fabian (1968) ^[10] demostraron que alcanzaríamos una tasa de convergencia óptima con (o ). Lai y Robbins ^{[11] [12]} diseñaron procedimientos adaptativos para estimar tal que tiene una varianza asintótica mínima. Sin embargo, la aplicación de tales métodos óptimos requiere mucha información a priori que es difícil de obtener en la mayoría de las situaciones. Para superar esta deficiencia, Polyak (1991) ^[13] y Ruppert (1988) ^[14] desarrollaron independientemente un nuevo algoritmo óptimo basado en la idea de promediar las trayectorias. Polyak y Juditsky ^[15] también presentaron un método para acelerar Robbins-Monro para problemas de búsqueda de raíces lineales y no lineales mediante el uso de pasos más largos y el promedio de las iteraciones. El algoritmo tendría la siguiente estructura: La convergencia de a la raíz única depende de la condición de que la secuencia de pasos disminuya lo suficientemente lentamente. Es decir ${\textstyle O(1/{\sqrt {n}})}$ ${\textstyle a_{n}=\bigtriangledown ^{2}f(\theta ^{*})^{-1}/n}$ ${\textstyle a_{n}={\frac {1}{(nM'(\theta ^{*}))}}}$ ${\textstyle M'(\theta ^{*})}$ ${\textstyle \theta _{n}}$ $${\displaystyle \theta _{n+1}-\theta _{n}=a_{n}(\alpha -N(\theta _{n})),\qquad {\bar {\theta }}_{n}={\frac {1}{n}}\sum _{i=0}^{n-1}\theta _{i}}$$ ${\displaystyle {\bar {\theta }}_{n}}$ ${\displaystyle \theta ^{*}}$ ${\displaystyle \{a_{n}\}}$

A1) $${\displaystyle a_{n}\rightarrow 0,\qquad {\frac {a_{n}-a_{n+1}}{a_{n}}}=o(a_{n})}$$

Por lo tanto, la secuencia con satisface esta restricción, pero no lo hace, de ahí los pasos más largos. Bajo los supuestos delineados en el algoritmo de Robbins-Monro, la modificación resultante dará como resultado la misma tasa de convergencia asintóticamente óptima pero con una política de tamaño de paso más robusta. ^[15] Antes de esto, la idea de usar pasos más largos y promediar las iteraciones ya había sido propuesta por Nemirovski y Yudin ^[16] para los casos de resolución del problema de optimización estocástica con objetivos convexos continuos y para problemas de punto de silla convexo-cóncavo. Se observó que estos algoritmos alcanzaban la tasa no asintótica . ${\textstyle a_{n}=n^{-\alpha }}$ ${\textstyle 0<\alpha <1}$ ${\textstyle \alpha =1}$ ${\textstyle O(1/{\sqrt {n}})}$ ${\textstyle O(1/{\sqrt {n}})}$

Un resultado más general se da en el Capítulo 11 de Kushner y Yin ^[17] al definir el tiempo interpolado , el proceso interpolado y el proceso normalizado interpolado como ${\textstyle t_{n}=\sum _{i=0}^{n-1}a_{i}}$ ${\textstyle \theta ^{n}(\cdot )}$ ${\textstyle U^{n}(\cdot )}$

$${\displaystyle \theta ^{n}(t)=\theta _{n+i},\quad U^{n}(t)=(\theta _{n+i}-\theta ^{*})/{\sqrt {a_{n+i}}}\quad {\mbox{for}}\quad t\in [t_{n+i}-t_{n},t_{n+i+1}-t_{n}),i\geq 0}$$Sea el promedio de iteración y el error normalizado asociado . ${\displaystyle \Theta _{n}={\frac {a_{n}}{t}}\sum _{i=n}^{n+t/a_{n}-1}\theta _{i}}$ ${\displaystyle {\hat {U}}^{n}(t)={\frac {\sqrt {a_{n}}}{t}}\sum _{i=n}^{n+t/a_{n}-1}(\theta _{i}-\theta ^{*})}$

Con el supuesto A1) y el siguiente A2)

A2) Existe una matriz de Hurwitz y una matriz simétrica y definida positiva tal que converge débilmente a , donde es la estatisolución de donde es un proceso de Wiener estándar. ${\textstyle A}$ ${\textstyle \Sigma }$ ${\textstyle \{U^{n}(\cdot )\}}$ ${\textstyle U(\cdot )}$ ${\textstyle U(\cdot )}$ $${\displaystyle dU=AU\,dt+\Sigma ^{1/2}\,dw}$$ ${\textstyle w(\cdot )}$

satisfecho y definir . Luego, para cada Error al analizar (SVG (MathML se puede habilitar a través del complemento del navegador): Respuesta no válida ("La extensión Math no se puede conectar a Restbase") del servidor "http://localhost:6011/en.wikipedia.org/v1/":): {\textstyle t} , ${\textstyle {\bar {V}}=(A^{-1})'\Sigma (A')^{-1}}$

$${\displaystyle {\hat {U}}^{n}(t){\stackrel {\mathcal {D}}{\longrightarrow }}{\mathcal {N}}(0,V_{t}),\quad {\text{where}}\quad V_{t}={\bar {V}}/t+O(1/t^{2}).}$$

El éxito de la idea del promedio se debe a la separación de la escala de tiempo de la secuencia original y la secuencia promediada , siendo la escala de tiempo de la primera más rápida. ${\textstyle \{\theta _{n}\}}$ ${\textstyle \{\Theta _{n}\}}$

Aplicación en optimización estocástica

Supongamos que queremos resolver el siguiente problema de optimización estocástica

$${\displaystyle g(\theta ^{*})=\min _{\theta \in \Theta }\operatorname {E} [Q(\theta ,X)],}$$donde es diferenciable y convexo, entonces este problema es equivalente a encontrar la raíz de . Aquí se puede interpretar como un cierto costo "observado" en función de los efectos aleatorios y elegidos . En la práctica, puede ser difícil obtener una forma analítica de , el método de Robbins-Monro logra generar una secuencia para aproximar si se puede generar , en la que la esperanza condicional de dado es exactamente , es decir, se simula a partir de una distribución condicional definida por ${\textstyle g(\theta )=\operatorname {E} [Q(\theta ,X)]}$ ${\displaystyle \theta ^{*}}$ ${\displaystyle \nabla g(\theta )=0}$ ${\displaystyle Q(\theta ,X)}$ ${\displaystyle \theta }$ ${\displaystyle X}$ ${\displaystyle \nabla g(\theta )}$ ${\displaystyle (\theta _{n})_{n\geq 0}}$ ${\displaystyle \theta ^{*}}$ ${\displaystyle (X_{n})_{n\geq 0}}$ ${\displaystyle X_{n}}$ ${\displaystyle \theta _{n}}$ ${\displaystyle \nabla g(\theta _{n})}$ ${\displaystyle X_{n}}$

$${\displaystyle \operatorname {E} [H(\theta ,X)|\theta =\theta _{n}]=\nabla g(\theta _{n}).}$$

Aquí hay un estimador imparcial de . Si depende de , en general no hay una forma natural de generar un resultado aleatorio que sea un estimador imparcial del gradiente. En algunos casos especiales cuando se aplican los métodos de IPA o de razón de verosimilitud, entonces se puede obtener un estimador imparcial del gradiente . Si se considera como un proceso aleatorio subyacente "fundamental" que se genera independientemente de , y bajo algunas condiciones de regularización para operaciones de intercambio de derivadas e integrales de modo que , entonces da la estimación imparcial del gradiente fundamental. Sin embargo, para algunas aplicaciones tenemos que usar métodos de diferencias finitas en los que tiene una expectativa condicional cercana a pero no exactamente igual a ella. ${\displaystyle H(\theta ,X)}$ ${\displaystyle \nabla g(\theta )}$ ${\displaystyle X}$ ${\displaystyle \theta }$ ${\displaystyle H(\theta ,X)}$ ${\displaystyle H(\theta ,X)}$ ${\displaystyle X}$ ${\displaystyle \theta }$ ${\displaystyle \operatorname {E} {\Big [}{\frac {\partial }{\partial \theta }}Q(\theta ,X){\Big ]}=\nabla g(\theta )}$ ${\displaystyle H(\theta ,X)={\frac {\partial }{\partial \theta }}Q(\theta ,X)}$ ${\displaystyle H(\theta ,X)}$ ${\displaystyle \nabla g(\theta )}$

Definimos entonces una recursión de forma análoga al método de Newton en el algoritmo determinista:

$${\displaystyle \theta _{n+1}=\theta _{n}-\varepsilon _{n}H(\theta _{n},X_{n+1}).}$$

Convergencia del algoritmo

El siguiente resultado proporciona condiciones suficientes para que el algoritmo converja: ^[18] ${\displaystyle \theta _{n}}$

C1) ${\displaystyle \varepsilon _{n}\geq 0,\forall \;n\geq 0.}$

C2) ${\displaystyle \sum _{n=0}^{\infty }\varepsilon _{n}=\infty }$

C3) ${\displaystyle \sum _{n=0}^{\infty }\varepsilon _{n}^{2}<\infty }$

C4) ${\displaystyle |X_{n}|\leq B,{\text{ for a fixed bound }}B.}$

C5) ${\displaystyle g(\theta ){\text{ is strictly convex, i.e.}}}$

$${\displaystyle \inf _{\delta \leq |\theta -\theta ^{*}|\leq 1/\delta }\langle \theta -\theta ^{*},\nabla g(\theta )\rangle >0,{\text{ for every }}0<\delta <1.}$$

Luego converge casi con seguridad. ${\displaystyle \theta _{n}}$ ${\displaystyle \theta ^{*}}$

A continuación se ofrecen algunas explicaciones intuitivas sobre estas condiciones. Supongamos que es una variable aleatoria uniformemente acotada. Si C2) no se cumple, es decir , entonces es una secuencia acotada, por lo que la iteración no puede converger a si la estimación inicial está demasiado alejada de . En cuanto a C3), observe que si converge a entonces ${\displaystyle H(\theta _{n},X_{n+1})}$ ${\displaystyle \sum _{n=0}^{\infty }\varepsilon _{n}<\infty }$ $${\displaystyle \theta _{n}-\theta _{0}=-\sum _{i=0}^{n-1}\varepsilon _{i}H(\theta _{i},X_{i+1})}$$ ${\displaystyle \theta ^{*}}$ ${\displaystyle \theta _{0}}$ ${\displaystyle \theta ^{*}}$ ${\displaystyle \theta _{n}}$ ${\displaystyle \theta ^{*}}$

$${\displaystyle \theta _{n+1}-\theta _{n}=-\varepsilon _{n}H(\theta _{n},X_{n+1})\rightarrow 0,{\text{ as }}n\rightarrow \infty .}$$Por lo tanto, debemos tener , y la condición C3) lo garantiza. Una elección natural sería . La condición C5) es una condición bastante estricta sobre la forma de ; proporciona la dirección de búsqueda del algoritmo. ${\displaystyle \varepsilon _{n}\downarrow 0}$ ${\displaystyle \varepsilon _{n}=1/n}$ ${\displaystyle g(\theta )}$

Ejemplo (donde es apropiado el método del gradiente estocástico)^[8]

Supongamos que , donde es diferenciable y es una variable aleatoria independiente de . Entonces depende de la media de , y el método del gradiente estocástico sería apropiado en este problema. Podemos elegir ${\displaystyle Q(\theta ,X)=f(\theta )+\theta ^{T}X}$ ${\displaystyle f}$ ${\displaystyle X\in \mathbb {R} ^{p}}$ ${\displaystyle \theta }$ ${\displaystyle g(\theta )=\operatorname {E} [Q(\theta ,X)]=f(\theta )+\theta ^{T}\operatorname {E} X}$ ${\displaystyle X}$ ${\displaystyle H(\theta ,X)={\frac {\partial }{\partial \theta }}Q(\theta ,X)={\frac {\partial }{\partial \theta }}f(\theta )+X.}$

Algoritmo de Kiefer-Wolfowitz

El algoritmo de Kiefer-Wolfowitz fue introducido en 1952 por Jacob Wolfowitz y Jack Kiefer ^[19] y fue motivado por la publicación del algoritmo de Robbins-Monro. Sin embargo, el algoritmo fue presentado como un método que estimaría estocásticamente el máximo de una función.

Sea una función que tiene un máximo en el punto . Se supone que es desconocida; sin embargo, ciertas observaciones , donde , se pueden realizar en cualquier punto . La estructura del algoritmo sigue un método de tipo gradiente, en el que las iteraciones se generan como ${\displaystyle M(x)}$ ${\displaystyle \theta }$ ${\displaystyle M(x)}$ ${\displaystyle N(x)}$ ${\displaystyle \operatorname {E} [N(x)]=M(x)}$ ${\displaystyle x}$

${\displaystyle x_{n+1}=x_{n}+a_{n}\cdot \left({\frac {N(x_{n}+c_{n})-N(x_{n}-c_{n})}{2c_{n}}}\right)}$

donde y son independientes. En cada paso, el gradiente de se aproxima de forma similar a un método de diferencia central con . Por lo tanto, la secuencia especifica la secuencia de anchos de diferencia finita utilizados para la aproximación del gradiente, mientras que la secuencia especifica una secuencia de tamaños de paso positivos tomados a lo largo de esa dirección. ${\displaystyle N(x_{n}+c_{n})}$ ${\displaystyle N(x_{n}-c_{n})}$ ${\displaystyle M(x)}$ ${\displaystyle h=2c_{n}}$ ${\displaystyle \{c_{n}\}}$ ${\displaystyle \{a_{n}\}}$

Kiefer y Wolfowitz demostraron que, si se cumplen ciertas condiciones de regularidad, entonces convergerá a con una probabilidad tal que , y posteriormente Blum ^[4] en 1954 demostró que converge a casi con seguridad, siempre que: ${\displaystyle M(x)}$ ${\displaystyle x_{n}}$ ${\displaystyle \theta }$ ${\displaystyle n\to \infty }$ ${\displaystyle x_{n}}$ ${\displaystyle \theta }$

• ${\displaystyle \operatorname {Var} (N(x))\leq S<\infty }$Para todos . ${\displaystyle x}$
• La función tiene un único punto máximo (mínimo) y es fuertemente cóncava (convexa) ${\displaystyle M(x)}$
  • El algoritmo se presentó por primera vez con el requisito de que la función mantenga una fuerte convexidad global (concavidad) sobre todo el espacio factible. Dado que esta condición es demasiado restrictiva para imponerla sobre todo el dominio, Kiefer y Wolfowitz propusieron que es suficiente imponer la condición sobre un conjunto compacto que se sabe que incluye la solución óptima. ${\displaystyle M(\cdot )}$ ${\displaystyle C_{0}\subset \mathbb {R} ^{d}}$
• La función satisface las condiciones de regularidad de la siguiente manera: ${\displaystyle M(x)}$
  • Existe y tal que ${\displaystyle \beta >0}$ ${\displaystyle B>0}$ $${\displaystyle |x'-\theta |+|x''-\theta |<\beta \quad \Longrightarrow \quad |M(x')-M(x'')|<B|x'-x''|}$$
  • Existe y tal que ${\displaystyle \rho >0}$ ${\displaystyle R>0}$ $${\displaystyle |x'-x''|<\rho \quad \Longrightarrow \quad |M(x')-M(x'')|<R}$$
  • Para cada , existe algún tal que ${\displaystyle \delta >0}$ ${\displaystyle \pi (\delta )>0}$ $${\displaystyle |z-\theta |>\delta \quad \Longrightarrow \quad \inf _{\delta /2>\varepsilon >0}{\frac {|M(z+\varepsilon )-M(z-\varepsilon )|}{\varepsilon }}>\pi (\delta )}$$
• Las secuencias seleccionadas deben ser secuencias infinitas de números positivos tales que ${\displaystyle \{a_{n}\}}$ ${\displaystyle \{c_{n}\}}$
  • ${\displaystyle \quad c_{n}\rightarrow 0\quad {\text{as}}\quad n\to \infty }$
  • ${\displaystyle \sum _{n=0}^{\infty }a_{n}=\infty }$
  • ${\displaystyle \sum _{n=0}^{\infty }a_{n}c_{n}<\infty }$
  • ${\displaystyle \sum _{n=0}^{\infty }a_{n}^{2}c_{n}^{-2}<\infty }$

Una elección adecuada de secuencias, según lo recomendado por Kiefer y Wolfowitz, sería y . ${\displaystyle a_{n}=1/n}$ ${\displaystyle c_{n}=n^{-1/3}}$

Desarrollos posteriores y cuestiones importantes

1. El algoritmo de Kiefer-Wolfowitz requiere que para cada cálculo de gradiente se simulen al menos valores de parámetros diferentes para cada iteración del algoritmo, donde es la dimensión del espacio de búsqueda. Esto significa que cuando es grande, el algoritmo de Kiefer-Wolfowitz requerirá un esfuerzo computacional sustancial por iteración, lo que conduce a una convergencia lenta. ${\displaystyle d+1}$ ${\displaystyle d}$ ${\displaystyle d}$
   1. Para solucionar este problema, Spall propuso el uso de perturbaciones simultáneas para estimar el gradiente. Este método requeriría solo dos simulaciones por iteración, independientemente de la dimensión . ^[20] ${\displaystyle d}$
2. En las condiciones requeridas para la convergencia, puede resultar difícil determinar la capacidad de especificar un conjunto compacto predeterminado que cumpla con una fuerte convexidad (o concavidad) y contenga la solución única. Con respecto a las aplicaciones del mundo real, si el dominio es bastante grande, estas suposiciones pueden ser bastante restrictivas y muy poco realistas.

Desarrollos futuros

Se ha desarrollado una extensa literatura teórica en torno a estos algoritmos, en relación con las condiciones de convergencia, las tasas de convergencia, las generalizaciones multivariadas y de otro tipo, la elección adecuada del tamaño de paso, los posibles modelos de ruido, etc. [ ^{21] [22]} Estos métodos también se aplican en la teoría de control , en cuyo caso la función desconocida que deseamos optimizar o de la que deseamos encontrar el cero puede variar en el tiempo. En este caso, el tamaño del paso no debe converger a cero, sino que debe elegirse de manera que siga la función. ^{[21] , 2.ª ed., capítulo 3} ${\displaystyle a_{n}}$

C. Johan Masreliez y R. Douglas Martin fueron los primeros en aplicar la aproximación estocástica a la estimación robusta . ^[23]

La principal herramienta para analizar algoritmos de aproximaciones estocásticas (incluidos los algoritmos Robbins-Monro y Kiefer-Wolfowitz) es un teorema de Aryeh Dvoretzky publicado en 1956. ^[24]

Véase también

• Descenso de gradiente estocástico
• Reducción de la varianza estocástica

Referencias

1. Toulis, Panos; Airoldi, Edoardo (2015). "Estrategias de estimación escalables basadas en aproximaciones estocásticas: resultados clásicos y nuevos conocimientos". Estadística y Computación . 25 (4): 781–795. doi :10.1007/s11222-015-9560-y. PMC  4484776 . PMID  26139959.
2. Le Ny, Jerome. "Introducción a los algoritmos de aproximación estocástica" (PDF) . Polytechnique Montreal . Teaching Notes . Consultado el 16 de noviembre de 2016 .
3. Robbins, H. ; Monro, S. (1951). "Un método de aproximación estocástica". Anales de estadística matemática . 22 (3): 400. doi : 10.1214/aoms/1177729586 .
4. Blum, Julius R. (1 de junio de 1954). "Métodos de aproximación que convergen con probabilidad uno". Anales de estadística matemática . 25 (2): 382–386. doi : 10.1214/aoms/1177728794 . ISSN  0003-4851.
5. Sacks, J. (1958). "Distribución asintótica de procedimientos de aproximación estocástica". Anales de estadística matemática . 29 (2): 373–405. doi : 10.1214/aoms/1177706619 . JSTOR  2237335.
6. Nemirovski, A. ; Juditsky, A.; Lan, G.; Shapiro, A. (2009). "Enfoque de aproximación estocástica robusta para la programación estocástica". Revista SIAM sobre optimización . 19 (4): 1574. doi :10.1137/070704277.
7. Complejidad del problema y eficiencia del método en optimización, A. Nemirovski y D. Yudin, Wiley -Intersci. Ser. Matemática discreta 15 John Wiley Nueva York (1983).
8. Introducción a la búsqueda y optimización estocástica: estimación, simulación y control, JC Spall, John Wiley Hoboken, NJ , (2003).
9. Chung, KL (1 de septiembre de 1954). "Sobre un método de aproximación estocástica". Anales de estadística matemática . 25 (3): 463–483. doi : 10.1214/aoms/1177728716 . ISSN  0003-4851.
10. Fabian, Vaclav (1968-08-01). "Sobre la normalidad asintótica en la aproximación estocástica". Anales de estadística matemática . 39 (4): 1327–1332. doi : 10.1214/aoms/1177698258 . ISSN  0003-4851.
11. Lai, TL; Robbins, Herbert (1979-11-01). "Diseño adaptativo y aproximación estocástica". Anales de estadística . 7 (6): 1196–1221. doi : 10.1214/aos/1176344840 . ISSN  0090-5364.
12. Lai, Tze Leung; Robbins, Herbert (1 de septiembre de 1981). "Consistencia y eficiencia asintótica de estimaciones de pendientes en esquemas de aproximación estocástica". Zeitschrift für Wahrscheinlichkeitstheorie und Verwandte Gebiete . 56 (3): 329–360. doi : 10.1007/BF00536178 . ISSN  0044-3719. S2CID  122109044.
13. Polyak, BT (1991). "Nuevos procedimientos de tipo aproximación estocástica. (En ruso.)". Automatización y control remoto . 7 (7).
14. Ruppert, David (1988). Estimadores eficientes a partir de un proceso Robbins-Monro de convergencia lenta (Informe técnico 781). Facultad de Investigación de Operaciones e Ingeniería Industrial de la Universidad de Cornell.
15. Polyak, BT; Juditsky, AB (1992). "Aceleración de la aproximación estocástica mediante promedio". Revista SIAM sobre control y optimización . 30 (4): 838. doi :10.1137/0330046.
16. Sobre la convergencia de Cezari del método de descenso más pronunciado para aproximar puntos de silla de funciones convexo-cóncavas, A. Nemirovski y D. Yudin, Dokl. Akad. Nauk SSR 2939 , (1978 (ruso)), Soviet Math. Dokl. 19 (1978 (inglés)).
17. Kushner, Harold; George Yin, G. (17 de julio de 2003). Aproximación estocástica y algoritmos recursivos | Harold Kushner | Springer. www.springer.com. ISBN 9780387008943. Recuperado el 16 de mayo de 2016 .
18. Bouleau, N. ; Lepingle, D. (1994). Métodos numéricos para procesos estocásticos. Nueva York: John Wiley. ISBN 9780471546412.
19. Kiefer, J.; Wolfowitz, J. (1952). "Estimación estocástica del máximo de una función de regresión". Anales de estadística matemática . 23 (3): 462. doi : 10.1214/aoms/1177729392 .
20. Spall, JC (2000). "Aproximación estocástica adaptativa mediante el método de perturbación simultánea". IEEE Transactions on Automatic Control . 45 (10): 1839–1853. doi :10.1109/TAC.2000.880982.
21. Kushner, HJ ; Yin, GG (1997). Algoritmos y aplicaciones de aproximación estocástica . doi :10.1007/978-1-4899-2696-8. ISBN 978-1-4899-2698-2.
22. Aproximación estocástica y estimación recursiva , Mikhail Borisovich Nevel'son y Rafail Zalmanovich Has'minskiĭ, traducido por Israel Program for Scientific Translations y B. Silver, Providence, RI: American Mathematical Society, 1973, 1976. ISBN 0-8218-1597-0 . 
23. Martin, R.; Masreliez, C. (1975). "Estimación robusta mediante aproximación estocástica". IEEE Transactions on Information Theory . 21 (3): 263. doi :10.1109/TIT.1975.1055386.
24. Dvoretzky, Aryeh (1956). "Sobre la aproximación estocástica". En Neyman, Jerzy (ed.). Actas del Tercer Simposio de Berkeley sobre Estadística Matemática y Probabilidad, 1954-1955, vol. I. University of California Press. págs. 39-55. MR  0084911.