Skip to main content

$$ \newcommand{\R}{\mathbb{R}} \newcommand{\N}{\mathbb{N}} \newcommand{\Z}{\mathbb{Z}} $$

Lecture 17: function approximation

Goal: architectures and algorithms for approximating optimal value and policy functions.

References

• Chapter 3: Bertsekas, D. P. & Tsitsiklis, J. N. Neuro-dynamic programming. Athena Scientific (1996).

Problem formulation

We’d like to solve MDP

\[ \min_u \text{E}[c(x,u)]\ \text{s.t.}\ x^+ \sim P(x,u) \]

with infinite state and action spaces: \(x\in X = \R^d\), \(u\in U = \R^m\). But in this case the value \(v : X\rightarrow\R\) and policy \(\pi:X\rightarrow\Delta(U)\) cannot be represented with a finite-dimensional vector.

Idea: to use a digital computer to approximate these functions, we need new architectures and algorithms.

Architectures for approximation

We seek to approximate \(f:\R^n\rightarrow\R^\ell\) using \(\tilde{f} : \R^n\times\Theta\rightarrow\R^\ell\) (the function \(f\) could be \(v^*\), \(\pi^*\), \(\pi\), or \(v^\pi\)) where \(\theta\in\Theta\) are parameters.

Q: what properties are desirable in the approximating function \(\tilde{f}\)?

Clearly want to provide a good approximation of \(f\), but ideally it would also be easy to find good values for the parameters \(\theta\).

Broadly, we can consider architectures that are linear or nonlinear in the parameters \(\theta\).

A linear architecture takes the form

\[ \tilde{f}(x;\theta) = \sum_{k=1}^K \theta_k b_k(x) \]

where \(\{b_k\}_{k=1}^K\) are termed basis functions (although they may not form a vector space basis).

For nonlinear architectures, artificial neural networks are popular. A single-layer “perceptron” takes the form

\[ \tilde{f}(x;\theta) = \sum_{k=1}^K \theta_k \sigma\left(\sum_{\ell=1}^L \theta_{k,\ell} z_\ell(x)\right) \]

where \(\sigma:\R\rightarrow\R\) is an activation function and \(z:X\rightarrow\R^L\) is termed a feature vector, e.g. a partial state observation or approximate sufficient statistic in a finite-dimensional vector space.

Two of the most popular activation functions are the logistic map and the rectified linear unit (ReLU) or hinge illustrated below.

from imports import *

x = np.linspace(-6, 6, 500)
logistic = 1 / (1 + np.exp(-x))
relu = np.maximum(0, x)

fig, axes = plt.subplots(1, 2, figsize=(9, 3.5))

# logistic
axes[0].plot(x, logistic, lw=2)
axes[0].axhline(0, color='k', lw=0.5, ls='--')
axes[0].axhline(1, color='k', lw=0.5, ls='--')
axes[0].set_yticks([0, 0.5, 1])
axes[0].set_xticks([-4, -2, 0, 2, 4])
axes[0].set_title(r'$\sigma(x) = \dfrac{1}{1+e^{-x}}$', fontsize=13)
axes[0].set_xlabel('$x$')
axes[0].set_ylabel(r'$\sigma(x)$')
axes[0].spines['left'].set_position('zero')
axes[0].spines['bottom'].set_position('zero')
axes[0].spines[['top', 'right']].set_visible(False)

# relu
axes[1].plot(x, relu, lw=2, color='tab:orange')
axes[1].set_yticks([0, 2, 4, 6])
axes[1].set_xticks([-4, -2, 0, 2, 4])
axes[1].set_title(r'$\sigma(x) = \max(0,\, x)$', fontsize=13)
axes[1].set_xlabel('$x$')
axes[1].set_ylabel(r'$\sigma(x)$')
axes[1].spines['left'].set_position('zero')
axes[1].spines['bottom'].set_position('zero')
axes[1].spines[['top', 'right']].set_visible(False)

plt.tight_layout()
show_with_alt(fig,
    "Two side-by-side plots of activation functions. "
    "Left: logistic (sigmoid) function sigma(x) = 1/(1+exp(-x)), "
    "an S-shaped curve ranging from 0 to 1, with dashed asymptotes at y=0 and y=1, "
    "passing through (0, 0.5). "
    "Right: ReLU function sigma(x) = max(0,x), "
    "zero for x<0 and linear with slope 1 for x>=0, forming an elbow at the origin."
)

An illustration of this 1-layer network is shown below.

from imports import *
import matplotlib.patches as mpatches

fig, ax = plt.subplots(figsize=(10, 5))
ax.set_aspect('equal')
ax.axis('off')

r_in, r_hid, r_out = 0.13, 0.27, 0.13
x_in, x_hid, x_out = 0.5, 2.6, 4.7
L, K = 3, 3
y_in  = np.linspace(1.0, -1.0, L)
y_hid = np.linspace(1.0, -1.0, K)
y_out = 0.0

# edges (drawn first so nodes sit on top)
for yi in y_in:
    for yk in y_hid:
        ax.plot([x_in, x_hid], [yi, yk], color='0.75', lw=1.0, zorder=1)
for yk in y_hid:
    ax.plot([x_hid, x_out], [yk, y_out], color='0.75', lw=1.0, zorder=1)

# nodes
def node(x, y, r, label='', fontsize=13, fc='white', ec='black'):
    ax.add_patch(mpatches.Circle((x, y), r, fc=fc, ec=ec, lw=1.8, zorder=3))
    if label:
        ax.text(x, y, label, ha='center', va='center', fontsize=fontsize, zorder=4)

for yi in y_in:
    node(x_in, yi, r_in, fc='black', ec='black')
for yk in y_hid:
    node(x_hid, yk, r_hid, r'$\sigma$', fontsize=14)
node(x_out, y_out, r_out, fc='black', ec='black')

# node labels
for i, yi in enumerate(y_in):
    ax.text(x_in - r_in - 0.12, yi, rf'$z_{i+1}(x)$',
            ha='right', va='center', fontsize=12)
ax.text(x_out + r_out + 0.12, y_out, r'$\tilde{f}(x;\theta)$',
        ha='left', va='center', fontsize=13)

# edge weight labels on clean (non-crossing) representative edges
# theta_{k,ell}: on the top horizontal edge z_1 -> hidden_1
mx_kl = (x_in + x_hid) / 2
ax.text(mx_kl, y_in[0] + 0.10, r'$\theta_{k,\ell}$',
        ha='center', va='bottom', fontsize=12,
        bbox=dict(boxstyle='round,pad=0.1', fc='white', ec='none', alpha=0.9))

# theta_k: on the top edge hidden_1 -> output
mx_k = (x_hid + x_out) / 2
my_k = (y_hid[0] + y_out) / 2
ax.text(mx_k, my_k + 0.10, r'$\theta_k$',
        ha='center', va='bottom', fontsize=12,
        bbox=dict(boxstyle='round,pad=0.1', fc='white', ec='none', alpha=0.9))

# column headings
for xp, lbl in [(x_in, 'feature\nvector'), (x_hid, 'hidden layer'), (x_out, 'output')]:
    ax.text(xp, 1.55, lbl, ha='center', va='bottom', fontsize=10, color='0.45')

ax.set_xlim(-1.4, 6.2)
ax.set_ylim(-1.8, 2.3)

plt.tight_layout()
show_with_alt(fig,
    "Diagram of a single-layer perceptron. "
    "Left column: three filled input nodes labeled z_1(x), z_2(x), z_3(x) (the feature vector). "
    "Center column: three hidden nodes each labeled sigma (sigmoid activation function), "
    "fully connected to all inputs with weights theta_{k,ell}. "
    "Right: one filled output node labeled f-tilde(x; theta), "
    "connected to all hidden nodes with weights theta_k. "
    "Column headings: 'feature vector', 'hidden layer', 'output'."
)

Fact: any continuous function \(f:\R^n\rightarrow\R^\ell\) can be approximated to arbitrary accuracy by a sufficiently large single-layer network, i.e. with sufficiently large \(K,L\) and the right choice of weights \(\theta\).

• Cybenko, G. Approximation by superpositions of a sigmoidal function. Math. Control Signals Systems 2, 303–314 (1989).

An amusing quote from the Neurodynamic Programming book in light of the explosion of deep architectures: > “in practice, the number of hidden layers is usually one or two, and almost never more than three”

Algorithms for approximation

Assuming there exists \(\theta^*\in\Theta\) for which \(\tilde{f}_{\theta^*} \approx f\), how can we find \(\theta^*\)?

Q: propose an algorithmic approach to approximate \(\theta^*\).

Consider the conceptual optimization problem \[ \min_{\theta\in\Theta} \| \tilde{f}_\theta - f \|. \] Since the function norm \(\| \cdot \|\) can’t be evaluated on a computer, we’ll assume given a finite dataset \(\{(x,f(x))\}_{x\in\Xi}\) and consider instead \[ \min_{\theta\in\Theta} \sum_{x\in\Xi} \| \tilde{f}_\theta(x) - f(x) \|^2 \] (the square is added to make the problem smooth). We can then apply standard algorithms to this finite-dimensional nonlinear program.

Example (linear case): Suppose

\[ \tilde{f}(x;\theta) = \sum_{k=1}^K \theta_k b_k(x) \]

so that \(D_\theta\tilde{f}(x;\theta) = [b_1(x),\dots,b_K(x)]\). Then

\[ \min_{\theta\in\Theta} \sum_{x\in\Xi} \| \tilde{f}_\theta(x) - f(x) \| \]

is a linear least-squares problem.

Q: solve this NLP.

From stationarity we have \[ \sum_x (\tilde{f}(x;\theta) - f(x)) D\tilde{f}(x;\theta) = \sum_x \left( \sum_{k=1}^K \theta_k b_k(x) - f(x) \right)[b_1(x),\dots,b_K(x)] \] which is affine in \(\theta\)! So we can rewrite the equation as \(A \theta + d = 0\) and solve \(\theta = -A^{-1} d\).

Example (nonlinear case): For a multilayer perceptron, the function is a composition

\[ \tilde{f}(x;\theta) = \sigma_N(W_N\sigma_{N-1}(W_{N-1}\cdots\sigma_1(W_1 x)\cdots)) \]

so the gradient with respect to \(\theta\) is a product

\[ D_\theta \tilde{f} = D\sigma_N W_N D\sigma_{N-1} W_{N-1}\cdots D\sigma_1 W_1. \]

To evaluate the derivative, it is efficient to take a forward pass to evaluate the derivatives and then a backward pass to compute the matrix-vector products – this is referred to as backpropagation. In modern toolkits the derivatives will be computed analytically using autodifferentiation.