Radial Basis Function

/assets/img/rbf-cool-ex.svg

Radial Basis Function (RBF) é uma função definida no espaço euclidiano onde o valor de cada depende apenas de sua distância até a origem. Por exemplo uma radial function \(\phi\) em duas dimensões tem a seguinte forma:

\[\phi ( x,y) \ =\ \varphi ( r) ,\ r=\sqrt{x^{2} +y^{2}}\]

A função é de base radial apenas se não houver variação quando houver rotação mantendo o ponto de origem fixado. Segue alguns exemplos de funções radiais:

\[\begin{matrix} f( x) =\sqrt{1+\sigma ^{2} +x^{2}} & & f( x) =e^{-\frac{1}{2\sigma ^{2}} \ x^{2}}\\ & & \\ f( x) =x^{\sigma ^{2}}\ln x & & f( x) =\frac{1}{1+\sigma ^{2} x^{2}} \end{matrix}\]

Comportamento da função

RBF com variação no \(\sigma\), \(\lambda\) e em \(c\)

\[f( x) =\lambda e^{-\frac{1}{2\sigma ^{2}} \ (x-c)^{2}}\]


Somatórios de RBF’s são utilizado para criar interpolação. Este processo de interpolação é muito similar ao de uma rede neural e em algumas literaturas pode-se encontrar RBF’s Neural networks.

Segue o exemplo abaixo que mostra de uma forma interativa como as funciona a somatória de RBF’s:

\[f( x) =\lambda e^{-\frac{1}{2\sigma ^{2}} \ (x-c)^{2}} + \lambda e^{-\frac{1}{2\sigma ^{2}} \ (x-c)^{2}}\]


Através de sua versatilidade podemos mostrar exemplos de como construir algumas variedades de função através das RBF`s:

  • Somatório de múltiplos RBF`s e sinais
\[f(x)=-e^{-\left(\frac{1}{0.8^{2}}( x-1)^{2}\right)} +e^{-\left(\frac{1}{0.8^{2}}( x-0)^{2}\right)} -e^{-\left(\frac{1}{0.8^{2}}( x+1)^{2}\right)}\]
/assets/img/rbf3g.svg
  • Somatório de RBF`s reproduzindo uma função linear
\[f(x)=+e^{-\left(\frac{1}{0.9^{2}}( x-1.2)^{2}\right)} -e^{-\left(\frac{1}{0.9^{2}}( x+1.2)^{2}\right)}\]
/assets/img/rbf-linear.svg
  • Somatório de RBF`s reproduzindo a função cosseno
\[f( x) =e^{-\left(\frac{1}{\sqrt{2}^{2}}( x+2\pi )^{2}\right)} -e^{-\left(\frac{1}{\sqrt{2}^{2}}( x+\pi )^{2}\right)} +e^{-\left(\frac{1}{\sqrt{2}^{2}}( x+0)^{2}\right)} -e^{-\left(\frac{1}{\sqrt{2}^{2}}( x-\pi )^{2}\right)} +e^{-\left(\frac{1}{\sqrt{2}^{2}}( x-2\pi )^{2}\right)}\]
/assets/img/rbf-cos.svg

Interpolação

Depois de entender um pouco melhor o comportamento de RBF`s podemos perceber que Acima forma mostrados alguns exemplos de como a RBF funcionam e como a combinação delas pode criar um aproximador de praticamente qualquer função. Porém ainda precisamos de uma forma de encontrar os parametros necessários da RBF, e para isso podemos montar um sistema. Segue o exemplo:

Supondo que temos a seguinte entrada e saída:

\[\begin{array}{ c c } x & f( x)\\ 2 & 3\\ 3 & 6\\ 4 & 5 \end{array}\]
/assets/img/rbfex1input.svg

E que para resolver essa regressão vamos utilizar 2 centros:

\[\begin{array}{l} C1=2\\ C2=4 \end{array}\]

Um sigma constante:

\[\sigma ^{2} =1\]

E a RBF:

\[R_{p} = e^{-\frac{1}{2\sigma ^{2}} .\parallel ( X_{i}) -( X_{p}) \parallel ^{2}}\]

Então podemos construir um sistema

\[\begin{cases} e^{-0.5\parallel ( x_{1}) -( c_{1}) \parallel ^{2}} a_{1} +e^{-0.5\parallel ( x_{1}) -( c_{2}) \parallel ^{2}} a_{2} +b=output_{1}\\ e^{-0.5\parallel ( x_{2}) -( c_{1}) \parallel ^{2}} a_{1} +e^{-0.5\parallel ( x_{2}) -( c_{2}) \parallel ^{2}} a_{2} +b=output_{2}\\ e^{-0.5\parallel ( x_{3}) -( c_{1}) \parallel ^{2}} a_{1} +e^{-0.5\parallel ( x_{3}) -( c_{2}) \parallel ^{2}} a_{2} +b=output_{3} \end{cases}\] \[\begin{cases} e^{-0.5\parallel ( 2) -( 2) \parallel ^{2}} a_{1} +e^{-0.5\parallel ( 2) -( 4) \parallel ^{2}} a_{2} +b=3\\ e^{-0.5\parallel ( 3) -( 2) \parallel ^{2}} a_{1} +e^{-0.5\parallel ( 3) -( 4) \parallel ^{2}} a_{2} +b=6\\ e^{-0.5\parallel ( 4) -( 2) \parallel ^{2}} a_{1} +e^{-0.5\parallel ( 4) -( 4) \parallel ^{2}} a_{2} +b=5 \end{cases}\\\] \[\\ \begin{cases} e^{-0.0} a_{1} +e^{-2.0} a_{2} +b=3\\ e^{-0.5} a_{1} +e^{-0.5} a_{2} +b=6\\ e^{-2.0} a_{1} +e^{-0.0} a_{2} +b=5 \end{cases}\\\] \[\\ \begin{cases} 1.0000a_{1} +0.1353 a_{2} +b1.0000=3\\ 0.6065a_{1} +0.6065 a_{2} +b1.0000=6\\ 0.1353a_{1} +1.0000 a_{2} +b1.0000=5 \end{cases}\]

Conseguimos resolver o sistema através do método da matriz inversa para matrizes quadradas e o da pseudo inversa para as demais.

\[\begin{array}{l} \underbrace{\begin{bmatrix} 1.0000 & 0.1353 & 1.0000\\ 0.6065 & 0.6065 & 1.0000\\ 0.1353 & 1.0000 & 1.0000 \end{bmatrix}}_{\begin{bmatrix} R \end{bmatrix}} .\underbrace{\begin{bmatrix} a1\\ a2\\ a3 \end{bmatrix}}_{\begin{bmatrix} a \end{bmatrix}} =\underbrace{\begin{bmatrix} 3\\ 6\\ 5 \end{bmatrix}}_{\begin{bmatrix} A \end{bmatrix}} \end{array}\] \[\\ \begin{bmatrix} R \end{bmatrix}_{3x3} .\begin{bmatrix} a \end{bmatrix}_{3x1} =\begin{bmatrix} A \end{bmatrix}_{3x1}\\\] \[\\ \begin{bmatrix} R \end{bmatrix} =\begin{bmatrix} 1.0000 & 0.1353 & 1.0000\\ 0.6065 & 0.6065 & 1.0000\\ 0.1353 & 1.0000 & 1.0000 \end{bmatrix}_{3x3}\\\] \[\\ \frac{1}{[ R]} =\left([ R]^{t} .[ R]\right)^{-1} .[ R]^{t}\\\] \[\\ [ R]^{t} .[ R] =\begin{bmatrix} 1.0000 & 0.6065 & 0.1353\\ 0.1353 & 0.6065 & 1.0000\\ 1.0000 & 1.0000 & 1.0000 \end{bmatrix}_{3x3} .\ \begin{bmatrix} 1.0000 & 0.1353 & 1.0000\\ 0.6065 & 0.6065 & 1.0000\\ 0.1353 & 1.0000 & 1.0000 \end{bmatrix}_{3x3}\\\] \[\\ [ R]^{t} .[ R] =\begin{bmatrix} 1.3862 & 0.6386 & 1.7419\\ 0.6386 & 1.3262 & 1.7419\\ 1.7419 & 1.7419 & 3.0000 \end{bmatrix}_{3x3}\\\] \[\\ \left([ R]^{t} .[ R]\right)^{-1} =\begin{bmatrix} 248.9582 & 247.6207 & -288.3246\\ 247.6207 & 248.9582 & -288.3246\\ -288.3246 & -288.3246 & 335.1485 \end{bmatrix}_{3x3}\\\] \[\\ \left([ R]^{t} .[ R]\right)^{-1} .[ R]^{t} \ =\begin{bmatrix} 248.9582 & 247.6207 & -288.3246\\ 247.6207 & 248.9582 & -288.3246\\ -288.3246 & -288.3246 & 335.1485 \end{bmatrix}_{3x3} .\ \begin{bmatrix} 1.0000 & 0.1353 & 1.0000\\ 0.6065 & 0.6065 & 1.0000\\ 0.1353 & 1.0000 & 1.0000 \end{bmatrix}_{3x3}\\\] \[\\ \left([ R]^{t} .[ R]\right)^{-1} .[ R]^{t} \ =\begin{bmatrix} -5.8546 & 12.8657 & -7.0111\\ -7.0111 & 12.8657 & -5.8546\\ 7.8034 & -146069 & 7.8034 \end{bmatrix}_{3x3} =\frac{1}{[ R]}\\\] \[\\ \begin{bmatrix} a \end{bmatrix} =\frac{1}{[ R]} \ \begin{bmatrix} A \end{bmatrix}\\ \\\] \[\begin{bmatrix} +24.5749\\ +26.8879\\ -25.2138 \end{bmatrix}_{3x1} =\begin{bmatrix} -5.8546 & 12.8657 & -7.0111\\ -7.0111 & 12.8657 & -5.8546\\ 7.8034 & -146069 & 7.8034 \end{bmatrix}_{3x3} .\begin{bmatrix} 3\\ 6\\ 5 \end{bmatrix}_{3x1}\\\] \[\\ f( x) =\sum ^{2}_{p=1} \lambda e^{-\frac{1}{2\sigma ^{2}} \ (x)^{2}}\\ \\\]

Substituímos os \(\lambda\) pelos valores encontrados e podemos testar nossa função de interpolação:

\[\\ f\ ( x) \ =a_{1} e^{-0.5\parallel ( X) -( 2) \parallel ^{2}} +a_{2} e^{-0.5\parallel ( X) -( 4) \parallel ^{2}} +b\\ \\ f\ ( x) \ =24.5749e^{-0.5\parallel ( X) -( 2) \parallel ^{2}} +26.8879e^{-0.5\parallel ( X) -( 4) \parallel ^{2}} -25.2138\\ \\ X=( 2) \Longrightarrow 24.5749e^{0.0} +26.8879e^{-2.0} -25.2138\ =3\\ \\ X=( 3) \Longrightarrow 24.5749e^{-0.5} +26.8879e^{-0.5} -25.2138\ =6\\ \\ X=( 4) \Longrightarrow 24.5749e^{-2.0} +26.8879e^{0.0} -25.2138\ =5\\\] \[\\ f( 1.5) =-2.3451\\ \\ f( 2.5) =+5.2027\\ \\ f( 3.5) =+6.4930\\ \\ f( 4.5) =-0.4055\\\]
/assets/img/rbfex1output.svg

Implementação

import numpy as np
from math import exp, pow
from numpy.linalg import inv
from functools import reduce
import matplotlib.pyplot as plt

def rbf(inp, out, center):
    def euclidean_norm(x1, x2):
        return sqrt(((x1 - x2)**2).sum(axis=0))

    def gaussian (x, c):
        return exp(-0.5 * pow(euclidean_norm(x, c), 2))

    R = np.ones((len(inp), (len(center) + 1)))

    for i, iv in enumerate(inp):
        for j, jv in enumerate(center):
            R[i, j] = (gaussian(inp[i], center[j]))

    Rt = R.transpose()
    RtR = Rt.dot(R)
    iRtR = inv(RtR)
    oneR = iRtR.dot(Rt)
    a = oneR.dot(out)

    def rbf_interpolation(x):
        sum = np.ones(len(center) + 1)

        for i, iv in enumerate(center):
            sum[i] = gaussian(x, iv)

        y = a * sum
        return reduce((lambda x, y: x + y), y)

    return rbf_interpolation

inp = np.array([2, 3, 4, 7])
out = np.array([3, 6, 5, 0])
center = np.array([0, 3, 7])

rbf_instance = rbf(inp, out, center)

input_test = input_test = np.linspace(0,10,100)
output_test = list(map(rbf_instance, input_test))

plt.plot(input_test, output_test)
plt.plot(inp, out, 'ro')
plt.ylabel('expected vs predicted')
plt.show()

/assets/img/rbf1.svg

inp = np.array([2, 3, 4, 5])
out = np.array([3, 1, 5, -2])
center = np.array([2, 3, 4])

rbf_instance = rbf(inp, out, center)

rbf_instance = rbf(inp, out, center)

input_test = input_test = np.linspace(2,5,100)
output_test = list(map(rbf_instance, input_test))

plt.plot(input_test, output_test)
plt.plot(inp, out, 'ro')
plt.ylabel('expected vs predicted')
plt.show()

/assets/img/rbf1.svg

inp = np.array([2, 4, 7, 9])
out = np.array([7, 8, 9, 5])
center = np.array([2, 3, 4])

rbf_instance = rbf(inp, out, center)

input_test = input_test = np.linspace(-5,15,100)
output_test = list(map(rbf_instance, input_test))

plt.plot(input_test, output_test)
plt.plot(inp, out, 'ro')
plt.ylabel('expected vs predicted')
plt.show()

/assets/img/rbf3.svg

XOR input

inp = np.array([np.array([1,1]), np.array([0,1]), np.array([0,0]), np.array([1,0])])
out = np.array([              0,               1,               0,             1])
center = np.array([ np.array([1,1]), np.array([0,0])])

rbf_instance = rbf(inp, out, center)


inp_test = np.array([np.array([1,1]), 
                     np.array([0,1]), 
                     np.array([0,0]), 
                     np.array([1,0])])
output = map(rbf_instance, inp_test)

def colorize(output):
    c = [None]* len(output)
    for i, iv in enumerate(output):
        if (output[i] > 0):
            c[i] = 'blue'
        else:
            c[i] = 'red'
    return c

inp_x = [1, 0, 0, 1]
inp_y = [1, 1, 0, 0]

c = colorize(output)

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

ax.scatter(inp_x, inp_y, output, color=c, depthshade=False)
plt.show()

/assets/img/rbf1.svg