**SVRAI** Métodos Aproximados Hasta ahora hemos supuesto que la función de valor tiene las características adecuadas para poder representarse como una tabla. Pero las tablas son representaciones útiles sólo para problemas pequeños y discretos. Los problemas con espacios de estado más grandes pueden requerir una cantidad inviable de memoria, y los métodos exactos que hemos visto en el tema anterior pueden llegar a requerir una cantidad excesivamente grande de recursos de computación (en tiempo y en espacio). En estos casos no queda más remedio que recurrir a métodos de programación dinámica aproximada, donde la solución puede no ser exacta. Una forma de aproximar las soluciones es utilizar la aproximación de la función de valor, que es en donde nos centraremos en este tema. Vamos a ofrecer varios enfoques para aproximar la función de valor y cómo incorporar la programación dinámica para derivar políticas aproximadamente óptimas. # Representaciones paramétricas Utilizaremos $U_\theta(s)$ para denotar una representación paramétrica de la función de valor, donde $\theta$ es el vector de parámetros. Hay muchas formas de representar $U_\theta(s)$, y veremos algunas de ellas un poco más adelante. Suponiendo que tenemos una aproximación de este tipo, podemos extraer una acción según: \begin{equation} \label{eq1} \pi(s) = \arg \max_a \Big( R(s, a) + \gamma \sum_{s'} T(s' | s, a)U_\theta(s')\Big) \end{equation} Las aproximaciones de la función de valor se utilizan a menudo en problemas con espacios de estado continuos, en cuyo caso la suma anterior puede sustituirse por una integral, y esta integral puede aproximarse utilizando muestras del modelo de transición. Una alternativa al cálculo de la ecuación eqn. [eq1] es aproximar la función de valor de la acción $Q(s, a)$. Si utilizamos $Q_\theta(s, a)$ para representar una aproximación paramétrica, podemos obtener una acción según: \begin{equation} \label{eq2} \pi(s) = \arg \max_a Q_\theta(s, a) \end{equation} En este tema vamos a ver cómo podemos aplicar la programación dinámica a un conjunto finito de estados $S = s_{1:m}$ para llegar a una aproximación paramétrica de la función de valor sobre el espacio de estados completo. Se pueden utilizar diferentes esquemas para generar este conjunto. Si el espacio de estados es de dimensión relativamente baja, podemos definir una cuadrícula/malla. Otro enfoque consiste en utilizar un muestreo aleatorio del espacio de estados. Sin embargo, algunos estados tienen más probabilidades de ser encontrados que otros y, por tanto, son más importantes a la hora de construir la función de valor. Podemos sesgar el muestreo hacia los estados más importantes realizando simulaciones con alguna política (quizás inicialmente aleatoria), a partir de un conjunto plausible de estados iniciales. Se puede utilizar un enfoque iterativo para mejorar nuestra aproximación de la función de valor en los estados de $S$. Alternamos entre la mejora de nuestras estimaciones de valor en $S$ mediante programación dinámica y el reajuste de nuestra aproximación en esos estados. El algoritmo $1$ proporciona una implementación en la que el paso de programación dinámica consiste en actualizaciones de Bellman tal y como se vio en la iteración de valores del tema anterior. Se puede crear un algoritmo similar para las aproximaciones del valor de la acción $Q_\theta$. !!!alg:Algoritmo 1. Iteración del valor aproximado para un MDP con la aproximación de la función de valor paramétrica `Uθ`. Realizamos copias de seguridad (como vimos en el tema anterior) en los estados de `S` para obtener un vector de utilidades `U`. A continuación, llamamos a `fit!(Uθ, S, U)` que, genéricamente, modifica la representación paramétrica `Uθ` para ajustar mejor el valor de los estados en `S` a las utilidades en `U`. Diferentes aproximaciones paramétricas tienen diferentes implementaciones para `fit!` ~~~~ struct ApproximateValueIteration Uθ # initial parameterized value function that supports fit! S # set of discrete states for performing backups k_max # maximum number of iterations end function solve(M::ApproximateValueIteration, 𝒫::MDP) Uθ, S, k_max = M.Uθ, M.S, M.k_max for k in 1:k_max U = [backup(𝒫, Uθ, s) for s in S] fit!(Uθ, S, U) end return ValueFunctionPolicy(𝒫, Uθ) end ~~~~ Todas las representaciones paramétricas que veremos en este tema siguen la estructura de este algoritmo, lo que puede cambiar es la evaluación de $U_\theta$ y la función que usemos en el ajuste de $U_\theta$ (`fit!`) a las estimaciones de las utilidades en los puntos de $S$. Podemos agrupar las representaciones paramétricas en dos categorías: 1. La primera categoría incluye los **métodos de aproximación local**, donde $\theta$ corresponde a los valores en los estados en $S$. Para evaluar $U_\theta(s)$ en un estado arbitrario $s$, tomamos una suma ponderada de los valores almacenados en $S$. 2. La segunda categoría incluye los **métodos de aproximación global**, donde $\theta$ no está directamente relacionada con los valores en los estados en $S$. De hecho, $\theta$ puede tener muchos menos, o incluso muchos más, componentes que los estados en $S$. Tanto la aproximación local como muchas aproximaciones globales pueden verse como una aproximación de función lineal $U_\theta(s) = \theta^⊤ \beta(s)$, donde los métodos difieren en cómo se define la función vectorial $\beta$. En los métodos de aproximación local, $\beta(s)$ determina cómo ponderar las utilidades de los estados en $S$ para aproximar la utilidad en el estado $s$. Las ponderaciones son generalmente no negativas y suman $1$. En muchos métodos de aproximación global, $\beta(s)$ se ve como un conjunto de funciones base que se combinan de forma lineal para obtener una aproximación para un $s$ arbitrario. También podemos aproximar la función de valor de la acción utilizando una función lineal, $Q_\theta(s, a) = \theta^⊤ \beta(s, a)$. En el contexto de las aproximaciones locales, podemos proporcionar aproximaciones sobre espacios de acción continuos eligiendo un conjunto finito de acciones $A \subseteq \mathscr{A}$. Nuestro vector de parámetros $\theta$ consistiría entonces en componentes $|S| \times |A|$, cada uno de los cuales correspondería a un par estado-acción. En este caso, la función $\beta(s, a)$ devolvería un vector con el mismo número de componentes que especifica cómo ponderar nuestro conjunto finito de valores estado-acción para obtener una estimación de la utilidad asociada al estado $s$ y a la acción $a$. # Vecinos más cercanos Un enfoque sencillo para la aproximación local en $s$ es utilizar el valor del estado de $S$ que es el **vecino más cercano** de $s$. Para utilizar este enfoque, necesitamos algo parecido a una métrica. Si usamos $d(s, s')$ para denotar la distancia entre dos estados $s$ y $s'$, entonces la función de valor aproximada es $U_\theta(s) = \theta_i$, donde $i = \arg \min_{j\in 1:m} d(s_j, s)$. Podemos generalizar este enfoque para promediar los valores de los **$k$ vecinos más cercanos**. Este enfoque sigue dando lugar a funciones de valor constantes a trozos, pero los diferentes valores de $k$ pueden dar lugar a mejores aproximaciones. El siguiente algoritmo proporciona una implementación de esta aproximación. !!!alg:Algoritmo 2. El método de $k$-vecinos más cercanos, que aproxima el valor de un estado `s` basándose en los `k` estados más cercanos de `S`, determinados por una función de distancia `d`. El vector `θ` contiene los valores de los estados de `S`. Se puede conseguir una mayor eficiencia utilizando estructuras de datos especializadas, como los _kd-trees_, implementados en `NearestNeighbors.jl`. ~~~~ mutable struct NearestNeighborValueFunction k # number of neighbors d # distance function d(s, s') S # set of discrete states θ # vector of values at states in S end function (Uθ::NearestNeighborValueFunction)(s) dists = [Uθ.d(s,s') for s' in Uθ.S] ind = sortperm(dists)[1:Uθ.k] return mean(Uθ.θ[i] for i in ind) end function fit!(Uθ::NearestNeighborValueFunction, S, U) Uθ.θ = U return Uθ end ~~~~ # Kernel de suavizado Otro método de aproximación local es por medio de un **kernel de suavizado**, en el que las utilidades de los estados en $S$ se suavizan sobre todo el espacio de estados. Este método requiere la definición de una función $k(s, s')$ que relaciona los pares de estados $s$ y $s'$. Por lo general, queremos que $k(s, s')$ sea mayor para los estados que están más cerca entre sí porque esos valores nos indican cómo ponderar juntas las utilidades asociadas a los estados en $S$. Este método da como resultado la siguiente aproximación lineal: \begin{equation} \label{eq3} U_\theta(s) = \sum_{i=1}^m \theta_i \beta_i(s) = \theta^⊤ \beta(s) \end{equation} donde \begin{equation} \label{eq4} \beta_i(s) = \frac{k(s, s_i)}{\sum_{j=1}^m k(s, s_j)} \end{equation} El siguiente algoritmo proporciona una implementación de este método. !!!alg:Algoritmo 3. Aproximación de la función de valor localmente ponderada definida por una función de núcleo `k` y un vector de utilidades `θ` en estados en `S`. ~~~~C mutable struct LocallyWeightedValueFunction k # kernel function k(s, s') S # set of discrete states θ # vector of values at states in S end function (Uθ::LocallyWeightedValueFunction)(s) w = normalize([Uθ.k(s,s') for s' in Uθ.S], 1) return Uθ.θ ⋅ w end function fit!(Uθ::LocallyWeightedValueFunction, S, U) Uθ.θ = U return Uθ end ~~~~ Hay muchas maneras de definir una función kernel. Podemos definir nuestro núcleo simplemente como la inversa de la distancia entre estados: \begin{equation} \label{eq5} k(s, s') = max(d(s, s'), \epsilon)^{-1} \end{equation} donde $\epsilon$ es una pequeña constante positiva para evitar dividir por $0$ cuando $s = s'$. La figura de la derecha muestra las aproximaciones de valores del mismo ejemplo anterior utilizando varias funciones de distancia. Como podemos ver, el suavizado del kernel puede dar lugar a aproximaciones suaves de la función de valor, en contraste con los $k$ vecinos más cercanos. Otro kernel común es el kernel gaussiano: \begin{equation} \label{eq6} k(s, s') = exp \Big(-\frac{d(s, s')^2}{2\sigma^2} \Big) \end{equation} donde $\sigma$ controla el grado de suavidad. # Interpolación lineal La **interpolación lineal** es otro enfoque común para la aproximación local. El caso unidimensional es sencillo, en el que el valor aproximado de un estado $s$ entre dos estados $s_1$ y $s_2$ es: \begin{equation} \label{eq7} U\theta(s) = \alpha \theta_1 + (1 -\alpha)\theta_2,\quad \quad \alpha = \frac{(s_2 - s)}{(s_2 - s_1)} \end{equation} La interpolación lineal puede extenderse a una malla multidimensional. En el caso bidimensional, llamado **interpolación bilineal**, interpolamos entre cuatro vértices. La interpolación bilineal se realiza a través de la interpolación lineal unidimensional, una vez en cada eje, siendo necesario conocer la utilidad de cuatro estados que ocupan los vértices de la malla. Dados cuatro vértices de coordenadas $s_1 = [x_1, y_1]$, $s_2 = [x_1, y_2]$, $s_3 = [x_2, y_1]$, y $s_4 = [x_2, y_2]$, y un estado $s = [x, y]$, el valor interpolado vendría dado por: \begin{equation} \begin{aligned} \label{eq8} U\theta(s) &= \alpha\theta_{12} + (1 -\alpha)\theta_{34}\\ &= \frac{x_2 - x}{x_2 - x_1} \theta_{12} + \frac{x - x_1}{x-2 - x_1} \theta_{34}\\ &= \frac{x_2 - x}{x_2-x_1} (\alpha\theta_1 + (1 -\alpha)\theta_2) + \frac{x - x_1}{x_2 - x_1} (\alpha\theta_3 + (1-\alpha)\theta_4)\\ &= \frac{x_2 - x}{x_2 - x_1}\Big(\frac{y_2 - y}{y_2-y_1} \theta_1 + \frac{y - y_1}{y_2 - y_1} \theta_2\Big) + \frac{x - x_1}{x_2 - x_1}\Big(\frac{y_2 - y}{y_2- y_1} \theta_3 +\frac{y - y_1}{y_2 - y_1} \theta_4\Big)\\ &= \frac{(x_2 - x)(y_2 - y)}{(x_2 - x_1)(y_2 - y_1)}\theta_1 + \frac{(x_2 - x)(y - y_1)}{(x_2 - x_1)(y_2 - y_1)}\theta_2 + \frac{(x - x_1)(y_2 - y)}{(x_2 - x_1)(y_2 - y_1)}\theta_3 + \frac{(x - x_1)(y - y_1)}{(x_2 - x_1)(y_2 - y_1)}\theta_4 \end{aligned} \end{equation} La interpolación resultante pondera cada vértice en función del área de su cuadrante opuesto, como muestran las siguientes figuras. La última de ellas muestra el resultado de aplicar este método al ejemplo visto anteriormente.    La **interpolación multilineal** en $d$ dimensiones se consigue de forma similar interpolando linealmente a lo largo de cada eje, lo que requiere $2^d$ vértices. También en este caso, la utilidad de cada vértice se pondera en función del volumen del hiperrectángulo opuesto. La interpolación multilineal se implementa en el algoritmo siguiente. !!!alg:Algoritmo 4. Un método para llevar a cabo la interpolación multilineal para estimar el valor del vector de estado, `s`, para valores de estado conocidos, `θ`, sobre una cuadrícula definida por un vértice inferior izquierdo, `o`, y un vector de anchuras, `δ`. Todos los vértices de la rejilla pueden escribirse `o + δ.*i` para algún vector de enteros no negativos `i`. El paquete `Interpolations.jl` de Julia también proporciona métodos de interpolación multilineal y otros. ~~~~ mutable struct MultilinearValueFunction o # position of lower-left corner δ # vector of widths θ # vector of values at states in S end function (Uθ::MultilinearValueFunction)(s) o, δ, θ = Uθ.o, Uθ.δ, Uθ.θ Δ = (s - o)./δ # Multidimensional index of lower-left cell i = min.(floor.(Int, Δ) .+ 1, size(θ) .- 1) vertex_index = similar(i) d = length(s) u = 0.0 for vertex in 0:2^d-1 weight = 1.0 for j in 1:d # Check whether jth bit is set if vertex & (1 << (j-1)) > 0 vertex_index[j] = i[j] + 1 weight *= Δ[j] - i[j] + 1 else vertex_index[j] = i[j] weight *= i[j] - Δ[j] end end u += θ[vertex_index...]*weight end return u end function fit!(Uθ::MultilinearValueFunction, S, U) Uθ.θ = U return Uθ end ~~~~ !!!note: Interpolación simplicial La interpolación multilineal puede ser ineficiente en dimensiones elevadas. En lugar de ponderar las contribuciones de $2^d$ puntos, la **interpolación simplicial** considera sólo $d + 1$ puntos en la vecindad de un estado dado para producir una superficie continua que coincida con los puntos de muestra conocidos. # Regresión Lineal Un enfoque de aproximación global simple es la **regresión lineal**, donde $U_\theta(s)$ es una combinación lineal de funciones base, también denominadas comúnmente **características**. Estas funciones base son generalmente una función no lineal del estado $s$ y se combinan en una función vectorial $\beta(s)$ (o $\beta(s, a)$), dando lugar a las aproximaciones: \begin{equation} \label{eq10} U_\theta(s) = \theta^⊤ \beta(s)\quad \quad Q_\theta(s, a) = \theta^⊤ \beta(s, a) \end{equation} Aunque nuestra aproximación es lineal con respecto a las funciones base, la aproximación resultante puede ser no lineal con respecto a las variables de estado subyacentes. La siguiente figura ilustra este concepto, donde se muestra un ejemplo en el que la regresión lineal con funciones base no lineales es lineal en dimensiones superiores. En este caso, la regresión polinómica puede considerarse lineal en un espacio tridimensional. La función existe en el plano formado por sus bases, pero no ocupa todo el plano porque los términos no son independientes.  Añadir más funciones base generalmente mejora la capacidad de igualar las utilidades objetivo en los estados en $S$, pero demasiadas funciones base pueden llevar a aproximaciones pobres en otros estados (lo que se conoce como **sobreajuste**). Existen metodologías para elegir un conjunto adecuado de funciones de base para nuestro modelo de regresión. El ajuste de los modelos lineales implica determinar el vector $\theta$ que minimiza el error cuadrático de las predicciones en los estados de $S = s_{1:m}$. Si las utilidades asociadas a esos estados se denotan como $u_{1:m}$, entonces queremos encontrar la $\theta$ que minimiza \begin{equation} \label{eq11} \sum_{i=1}^m ( \hat{U}_\theta(s_i) - u_i)^2 = \sum_{i=1}^m (\theta^⊤ \beta(s_i) - u_i)^2 \end{equation} La $\theta$ óptima puede calcularse mediante algunas operaciones matriciales sencillas. !!!note:Cálculo de $\theta$ Primero construimos una matriz $\mathbf{X}$ donde cada una de las $m$ filas $\mathbf{X_{i,:}}$ contiene $\beta(s_i)^⊤$. Se puede demostrar que el valor de $\theta$ que minimiza el error cuadrático es \begin{equation} \label{eq12} \theta = \Big(\mathbf{X}^⊤ \mathbf{X}\Big)^{-1} \mathbf{X}^⊤ u_{1:m} = \mathbf{X}^+ u_{1:m} \end{equation} donde $\mathbf{X}^+$ es la pseudoinversa de Moore-Penrose de la matriz $\mathbf{X}$. La pseudoinversa suele implementarse calculando primero la descomposición de valores singulares, $\mathbf{X} = \mathbf{U}\mathbf{\Sigma}\mathbf{U}^*$, de donde tenemos \begin{equation} \label{eq13} \mathbf{X}^+ = \mathbf{U}\mathbf{\Sigma}^+\mathbf{U}^* \end{equation} La pseudoinversa de la matriz diagonal $\mathbf{\Sigma}$ se obtiene tomando el recíproco de cada elemento no nulo de la diagonal y luego transponiendo el resultado. La siguiente figura muestra cómo se ajustan las utilidades de los estados en $S$ con varias familias de funciones base. Diferentes elecciones de funciones base dan lugar a diferentes errores.  El siguiente algoritmo proporciona una implementación para evaluar y ajustar modelos de regresión lineal de la función de valor. !!!alg:Algoritmo 6. Aproximación de la función de valor de regresión lineal, definida por una función de vector base `β` y un vector de parámetros `θ`. La función `pinv` implementa la psuedoinversa. Julia y otros lenguajes soportan la operación de barra invertida, que nos permite escribir `X \ U` en lugar de `pinv(X)*U` en la función `fit!`. ~~~~ mutable struct LinearRegressionValueFunction β # basis vector function θ # vector of parameters end function (Uθ::LinearRegressionValueFunction)(s) return Uθ.β(s) ⋅ Uθ.θ end function fit!(Uθ::LinearRegressionValueFunction, S, U) X = hcat([Uθ.β(s) for s in S]...)' Uθ.θ = X / U # pinv(X)*U return Uθ end ~~~~ # Regresión por redes neuronales La **regresión por redes neuronales** nos libera de tener que construir un conjunto adecuado de funciones base, como se requiere en la regresión lineal. En su lugar, se utiliza una red neuronal para representar nuestra función de valor. La entrada de la red neuronal serían las variables de estado, y la salida sería la estimación de la utilidad. Los parámetros $\theta$ corresponderían a los pesos de la red neuronal. Podemos optimizar los pesos de la red para conseguir un objetivo concreto. En el contexto de la programación dinámica aproximada, querríamos minimizar el error de nuestras predicciones, tal como hicimos en la sección anterior. Sin embargo, la minimización del error cuadrático no puede hacerse mediante simples operaciones matriciales. En su lugar, generalmente tenemos que confiar en técnicas de optimización como el **descenso del gradiente**. Afortunadamente, el cálculo del gradiente de las redes neuronales puede hacerse exactamente mediante la aplicación directa de la regla de la cadena de derivadas. # Resumen * Para problemas grandes o continuos, podemos intentar encontrar políticas aproximadas representadas por modelos paramétricos de la función de valor. * Los enfoques adoptados en este tema implican la aplicación iterativa de pasos de programación dinámica en un conjunto finito de estados y el refinamiento de nuestra aproximación paramétrica. * Las técnicas de aproximación local aproximan la función de valor basándose en los valores de los estados cercanos con valores conocidos. * Una variedad de técnicas de aproximación local incluyen el vecino más cercano, el suavizado del kernel, la interpolación lineal y la interpolación simplicial. * Las técnicas de aproximación global incluyen la regresión lineal y la regresión por redes neuronales. * Se pueden obtener funciones de utilidad no lineales cuando se utiliza la regresión lineal, si se combina con una selección adecuada de funciones base no lineales. * La regresión por redes neuronales nos libera de tener que especificar las funciones base, pero su ajuste es más complejo y generalmente requiere que utilicemos el descenso de gradiente para ajustar nuestra aproximación paramétrica de la función de valor. # Ejercicios **Ejercicio 1.** Dados los estados $s_1 = (4, 5)$, $s_2 = (2, 6)$, y $s_3 = (-1, -1)$, y sus correspondientes valores, $U(s_1) = 2$, $U(s_2) = 10$, y $U(s_3) = 30$, calcula el valor en el estado $s = (1, 2)$ utilizando la aproximación local del vecino más cercano con una métrica de distancia $L_1$, con una métrica de distancia $L_2$ y con una métrica de distancia $L_\infty$. **Ejercicio 2.** Queremos estimar el valor en un estado $s$ dados los valores en un conjunto de dos estados $S = \{s_1, s_2\}$. Si queremos utilizar la iteración de valores por aproximación local, ¿cuáles de las siguientes funciones de ponderación son válidas? Si no son válidas, ¿cómo se podrían modificar las funciones de ponderación para hacerlas válidas? * $\beta(s) = [1, 1]$ * $\beta(s) = [1 - \lambda, \lambda]$ donde $\lambda \in [0, 1]$ * $\beta(s) = [e^{(s-s_1)^2}, e^{(s-s_2)^2}]$