Modelos de Vetores de Suporte

Uma Introdução ao Aprendizado Estatístico de Máquina

SUMÁRIO

🗒️ Resumo

🗒️ 1. Introdução

🗒️ 2. Fundamentação teórica

🗒️ 3. Primeiros passos Hands-on

🗒️ 4. Aplicações com SVM

• 📋 4.1 Linguagem R

• 📋 4.2 Linguagem Python

🗒️ 5. Outros tipos de Modelos de Vetores de Suporte

🗒️ 6. Exercícios

🗒️ Apêndice

5. Algumas Extensões dos Modelos de Vetores de Suporte

5.1 Classificação

5.1.1. Máquina de Vetores de Suporte probabilística

O SVM tradicional calcula uma função de decisão \(f(\textbf{x})\) tal que \(\text{sign}(f(\textbf{x}))\) pode ser utilizado para prever a classe de uma nova observação. Porém, muitas aplicações necessitam de uma probabilidade posteriori (probabilidade revisada de uma ocorrência após a consideração de novas evidências ou informações) \(P(Y = 1|\textbf{x})\) ao invés de um rótulo (Lin, Lin, e Weng 2007). Platt et al. (1999) propôs uma aproximação para a probabilidade posteriori dada por uma função sigmoide

\[ \operatorname{P}(Y=1 \mid \textbf{x}) \approx P_{A, B}(f) \equiv \frac{1}{1+\exp (A f+B)}, ~~~\text { em que } f=\mathrm{f}(\textbf{x}). \]

A melhor configuração de parâmetros \(z^*=\left(A^{*}, B^{*}\right)\) é determinada resolvendo o seguinte problema via Máxima Verossimilhança (com \(N_+\) a quantidade dos \(y_{i}\) positivos e \(N_-\) a quantidade dos \(y_{i}\) negativos) (Lin, Lin, e Weng 2007).

Os parâmetros \(A\) e \(B\) são ajustados utilizando a estimador de de máxima verossimilhança de um conjunto de treinamento \((f_i, y_i)\). Os parâmetros A e B são encontrados minimizando a probabilidade de log negativo da verossimilhança dos dados de treinamento, que pode ser visto como uma função de erro da entropia cruzada:

\[ \min _{z=(A, B)} F(z)=-\sum_{i=1}^{n}\left(t_{i} \log \left(p_{i}\right)+\left(1-t_{i}\right) \log \left(1-p_{i}\right)\right) \text {, } \] para \(p_{i}=P_{A, B}\left(f_{i}\right)\), e \(t_{i}=\left\{\begin{array}{ll}\frac{N_{+}+1}{N_{+}+2} & \text { if } y_{i}=+1 \\ \frac{1}{N_{-}+2} & \text { if } y_{i}=-1\end{array}, i=1, \ldots, n\right.\).

A minimização é no espaço dos dois parâmetros e portanto, pode ser alcanzada usando qualquer algoritmo de otimização, como por exemplo, algoritmos tipo gradiente como o algoritmo de Newton (Nocedal e Wright 1999).

Os pacotes citados ao longo destes trabalho para R e Python já possuem as implementações desta metodologia. No kernlab, por exemplo, basta adicionar o argumento prob.model=T na função ksvm. Similarmente, na função SVC do sklearn.svm do Python, basta adicionar o argumento probability=True. Abaixo a metodologia é explorada na Linguagem R.

require(kernlab)
data(spam)
## estimação do modelo
modelo <- ksvm(type~.,data=spam, kernel="rbfdot",
               kpar=list(sigma=0.05),C=5, prob.model=T)

modelo

## Support Vector Machine object of class "ksvm" 
## 
## SV type: C-svc  (classification) 
##  parameter : cost C = 5 
## 
## Gaussian Radial Basis kernel function. 
##  Hyperparameter : sigma =  0.05 
## 
## Number of Support Vectors : 1547 
## 
## Objective Function Value : -2082.468 
## Training error : 0.022169 
## Probability model included.

## predição para estimativa probabilística
preds <- predict(modelo,spam[,-58],type="prob")
head(preds)

##           nonspam      spam
## [1,] 0.0232441366 0.9767559
## [2,] 0.0001050739 0.9998949
## [3,] 0.0502387679 0.9497612
## [4,] 0.0169331336 0.9830669
## [5,] 0.0170070945 0.9829929
## [6,] 0.0493621241 0.9506379

5.1.2 Máquina de Vetores de Suporte para mais de duas classes

Quando a variável resposta assume três categorias ou mais, ajusta-se então um maior número de SVMs.

Existem dois métodos comuns (James et al. 2013): i) um-versus-um e ii) um-versus-todos. A Figura 1 abaixo exibe de forma esquemaática o problema de classificação para três categorias.

Figura 1. Classificação multiclasse. A: três categorias; B: método um-versus-um; C: método um-versus-todos. Fonte: Ng, Tse, e Tsui (2014)

A abordagem um-versus-um, considerando \(k\) classes (nosso caso \(k=3\)), iria exigir \({k \choose 2}\) ajustes de modelos diferentes, sendo uma abordagem bastante precisa e popularmente utilizada.

Porém, a abordagem um-versus-todos é menos onerosa computacionalmente uma vez que exige apenas \(k\) ajustes de modelos. Desta forma, para o método um-versus-todos, cada uma das classes é comparada as classes remanescentes. Sejam \(\mathbf{w}_k=(w_{0k},{w}_{1k},\ldots,{w}_{pk})\) os pesos resultantes de um ajuste SVM comparando a \(k\)-ésima classe (codificada como \(+1\)) com as outras (codificadas como \(-1\)) e \(\textbf{x}^{*}=\left(1,x_{1}^{*},\dots,x_{p}^{*}\right)^{'}\) uma observação-teste. Atribui-se a observação à classe para a qual \(w_{0k} + w_{1k}x^{*}_{1}+\cdots+w_{pk}x^{*}_{p}\) é maior, pois isso equivale a um alto nível de confiança que a observação de teste pertence à \(k\)-ésima classe e não a qualquer outra das classes.

As funções kernlab e o sklearn em R e Python respectivamnete implementam automaticamente SVM multiclasse por meio de uma estratégia um-versus-um e combinam vários SVMs binários.

5.2 Modelo de Vetores de Suporte para Regressão

Supondo um contexto no qual a variável resposta \(Y\) é contínua, \(y \in \mathbb{R}\) e as observações $ ^p $ é o vetor de variáveis explicativas, os modelos de vetores de suporte para regressão (SVR - Support Vector Regression) (Drucker et al. 1997) foram propostas como uma versão generalizada de máquinas de vetores de suporte. Da mesma forma, o método herda propriedades importantes da versão de classificação, uma vez que está baseado em vetores de suporte e no truque do kernel para lidar com cenários onde os dados não possuem comportamento linear.

Consideremos inicialmente a função linear:

\[ \tag{1}\label{svr:eq1} f(\textbf{x}_i) = \textbf{w} \cdot \textbf{x} + b \]

em que \(f(\textbf{x}_i)\) é determinada a partir da minimização da norma:

\[ \tag{2}\label{svr:eq2} \frac{\parallel \textbf{w} \parallel^{2}}{2} = \frac{\textbf{w} \cdot \textbf{w} }{2} \]

sujeita às condições:

\[ |y_{i} - f(x_{i})| \leq \epsilon, \; \forall \; i = 1, \dots, n \]

Tais condições advém do uso da função perda \(\epsilon\)-insensível (V. N. Vapnik 2000):

\[ \tag{3}\label{svr:eq3} L_{\epsilon}(f,y)= \begin{cases} 0, & \text{ se } | f(x)-y)|< \epsilon \\ | f(x)-y)| - \epsilon & \text{ caso contrário } \end{cases} \]

Então a ideia se estende do SVM para classificação, porém o que muda é a ótica de análise, pois desta vez, o interesse é ver como os dados se comportam dentro de um hipertubo (um hipertubo é uma margem que separa duas classes de dados num espaço multidimensional). A Figura 2, ilustra como os dados se comportariam num SVR com margens rígidas, observe que desta vez todos os dados se localizam dentro de um hipertubo. A Equação \(\eqref{svr:eq3}\) é conhecida como \(\epsilon\text{-insensível}\), e significa que o problema de otimização da Equação \(\eqref{svr:eq2}\) se restringe aos resíduos absolutos inferiores a \(\epsilon\), ou seja que alguns pontos podem não atender esta restrição.

Figura 2. SVR de margens rigidas.

Entretanto, como o interesse é obter um método mais flexível, de maneira análoga ao SVM com margens flexíveis para classificação, é possível a introdução de variáveis de folga \(\xi_{i} ,\xi^{\ast}_{i}\) para cada ponto. Sendo a ideia agora minimizar

\[\begin{equation}\tag{4}\label{svr:eq4} \frac{\textbf{w} \cdot \textbf{w} }{2} + C \sum^{n}_{i=1}(\xi_{i} +\xi^{\ast}_{i} ) \end{equation}\]

sujeita às seguintes restrições

\[\begin{eqnarray} y_{i} - f(x_{i}) \leq \epsilon + \xi_{i} \tag{5}\label{srv:eq5}\\ f(x_{i})- y_{i} \leq \epsilon + \xi^{\ast}_{i} \tag{6}\label{srv:eq6} \end{eqnarray}\] sendo \(\xi_{i} ,\xi^{\ast}_{i} \geq 0, \; \forall \; i = 1, \dots, n\), e \(C > 0\) é uma constante que permite controlar a penalidade imposta às observações que se encontram fora da margem \(\epsilon\). A Figura 3 ilustra o comportamento de dados dentro de um SVR com margens suaves, desta vez porém é possível que alguns pontos saiam do hipertubo, uma das intenções desta flexibilidade é evitar que aconteça o overfitting.

Figura 3. SVR de margens suaves.

A função de perda \(L_{\epsilon}(f,y)\) ignora os erros inferiores a \(\epsilon\) , tratando-os como nulos. A Figura 4 ilustra o comportamento dessa função de perda.

Figura 4. Função de perda \(\epsilon\)-insensível.

Com base na função \(\textbf{w} \cdot \textbf{w}/2\) e suas restrições, o Lagrangiano obtido a partir da expressões \(\eqref{srv:eq5}\) e \(\eqref{srv:eq6}\) utilizando multiplicadores de Lagrange em \(\boldsymbol{\alpha}, \boldsymbol{\alpha}^{\ast},\boldsymbol{\beta}, \boldsymbol{\beta}^{\ast}\) é dado por (Hamel 2011):

\[\begin{eqnarray} \max_{\boldsymbol{\alpha} \boldsymbol{\alpha}^{\ast}\boldsymbol{\beta} \boldsymbol{\beta}^{\ast}} \min_{ b \textbf{w},\boldsymbol{\xi} \boldsymbol{\xi}^{\ast}} L(\boldsymbol{\alpha}, \boldsymbol{\alpha}^{\ast},\boldsymbol{\beta}, \boldsymbol{\beta}^{\ast}, b ,\textbf{w},\boldsymbol{\xi} ,\boldsymbol{\xi}^{\ast}) = \max_{\boldsymbol{\alpha} \boldsymbol{\alpha}^{\ast},\boldsymbol{\beta} \boldsymbol{\beta}^{\ast}} \min_{ b ,\textbf{w},\boldsymbol{\xi} ,{\xi}^{\ast}} \Bigg( \frac{\textbf{w} \cdot \textbf{w} }{2}+ \nonumber \\ + C \sum^{n}_{i=1}(\xi_{i} +\xi^{\ast}_{i} ) - \sum^{n}_{i=1} \alpha_{i} (\xi_{i} + \epsilon - y_{i}+ \widehat{f}(\textbf{x}_{i})) + \\ - \sum^{n}_{i=1} \alpha^{\ast}_{i} (\xi^{*}_{i} + \epsilon - \widehat{f}(\textbf{x}_{i})+ y_{i}) - \nonumber \\ +\sum^{n}_{i=1}\beta_{i} \xi_{i}-\sum^{n}_{i=1}\beta^{\ast}_{i} \xi_{i}^{\ast} \Bigg), \end{eqnarray}\] em que \(\alpha_{i}, \alpha_{i}^{\ast},\beta_{i}, \beta_{i}^{\ast} \geq 0\), \(\widehat{f}(\textbf{x}_{i})= \textbf{w} \cdot \textbf{x} - b\).

Então, para encontrar os valores de \(\widehat{f}(\textbf{x})= \textbf{w} \cdot \textbf{x} - b\) é utilizada a seguinte versão primal do problema otimização:

\[ \tag{7}\label{svr:eq7} \max_{\boldsymbol{\alpha}^{\ast}, \boldsymbol{\alpha}} \Bigg[- \frac{1}{2} \sum^{n}_{i=1} \sum^{n}_{j=1}(\alpha^{\ast}_{i}-\alpha_{i})(\alpha^{\ast}_{j}-\alpha_{j})x_{i}x_{j} +\sum^{n}_{i=1} \alpha_{i}(y_{i} - \xi) - \alpha^{\ast}_{i}(y_{i} + \xi) \Bigg] \] em que \(0 \leq \alpha_{i} , \alpha^{\ast}_{i}\leq C, \: \forall i, \ldots, n\).

As derivadas de \(L=L(\boldsymbol{\alpha}, \boldsymbol{\alpha}^{\ast},\boldsymbol{\beta}, \boldsymbol{\beta}^{\ast}, b ,\textbf{w},\boldsymbol{\xi} ,\boldsymbol{\xi}^{\ast})\) são:

\[\begin{eqnarray} \frac{\partial L}{\partial b} &= \sum^{n}_{i=1}(\alpha^{\ast}_{i}-\alpha_{i}) &= 0 ,\tag{8}\label{svr:eq8} \\ \frac{\partial L}{\partial \textbf{w}} &= \textbf{w} - \sum^{n}_{i=1}(\alpha^{\ast}_{i}-\alpha_{i})x_{i} &= 0 ,\tag{9}\label{svr:eq9} \\ \frac{\partial L}{\partial \xi_{i}} &= C-\alpha_{i}-\beta_{i} &= 0, \tag{10}\label{svr:eq10}\\ \frac{\partial L}{\partial \xi^{\ast}_{i}} &= C-\alpha^{\ast}_{i}-\beta_{i} &= 0 \tag{11}\label{svr:eq11} \end{eqnarray}\]

Utilizando os resultados das Equações \(\eqref{svr:eq8}\)), \(\eqref{svr:eq9}\), \(\eqref{svr:eq10}\), \(\eqref{svr:eq11}\) na otimização \(\ref{svr:eq7}\). Temos que \[ \tag{12}\label{svr:eq12} \textbf{w}^{\ast}= \sum^{n}_{i=1}(\alpha^{\ast}_{i}-\alpha_{i})x_{i} \]

Substituindo \(\eqref{svr:eq12}\) em \(\eqref{svr:eq1}\) temos que:

\[ \tag{13}\label{svr:eq13} f(x_{i})= b + \sum^{n}_{i=1}(\alpha^{\ast}_{i}-\alpha_{i})x_{i} \cdot \textbf{x} \]

em que o valor de \(b\) é dado por:

\[ b^{\ast} = \frac{\sum^{n}_{i=1} \textbf{w} \cdot \textbf{x} -y_{i}}{n} \]

\[ \tag{14}\label{svr:eq14s} b^{\ast} = \frac{\sum^{n}_{i=1} \textbf{w} \cdot \textbf{w} -y_{i}}{n} \]

Para o caso não linear tem-se que:

\[ f(\textbf{x}_{i})= \textbf{w} \cdot \phi (\textbf{x}_{i})+b \] sendo que \(\phi\) é a função de mapeamento que mapeia a entrada \(\textbf{x}\) para o espaço \(\mathcal{R}^{n}\), para o espaço característico desejado. A Figura 5, ilustra como seria a transição dos dados para um espaço de característica superior.

Figura 5. Representação do hipertubo no espaço de características.

Colocando \(\phi\) na Equação \(\eqref{svr:eq14s}\), temos:

\[\begin{eqnarray}%\label{svr:eq15} f(x_{i})&=& b + \sum^{n}_{i=1}(\alpha^{\ast}_{i}-\alpha_{i}) \phi( x_{i} ) \cdot \phi (\textbf{x}) \label{svr:eq15} \\ &=& b + \sum^{n}_{i=1}(\alpha^{\ast}_{i}-\alpha_{i})K(x_{i},\textbf{x}) \end{eqnarray}\]

Por analogia aos resultados anteriores, o seguinte problema remete ao procedimento de otimização:

\[\begin{eqnarray} \max_{\alpha,\alpha^{\ast}} -\frac{1}{2} \sum^{n}_{i=1} \sum^{n}_{j=1}(\alpha^{\ast}_{i}-\alpha_{i})(\alpha^{\ast}_{j}-\alpha_{j})K(x_{i}x_{j}) - \epsilon \sum^{n}_{i=1}(\alpha^{\ast}_{i}-\alpha_{i})+ \sum^{n}_{i=1}y_{i}(\alpha^{\ast}_{i}-\alpha_{i}) \\ \text{sujeito a } \sum^{n}_{i=1}(\alpha^{\ast}_{i}-\alpha_{i})=0, 0\leq \alpha_{i} \leq C, \\ 0\leq \alpha_{i}^{\ast} \leq C \nonumber \end{eqnarray}\]

Ressaltando que os valores de \(\epsilon\) e \(C\) podem ser encontrados via processo de tunagem.

Um rápido exemplo de como utilizar o modelo de máquina de vetores de suporte para regressão pode ser visto abaixo em linguagem R. Para a linguagem Python o ajuste é realizado através do módulo sklearn.svm.SVR.

### PACOTES UTILIZADOS ---------------------------------------
require(kernlab) 

### GERANDO DATASET ------------------------------------------
set.seed(100)
n=1000
x=seq(-5,5,length.out = n)
y=x+sin(x)+rnorm(n,0,0.5)

plot(x,y,col='gray')

### COMPARAÇÃO REGRESSÃO C=1 e C=5 via holdout repetido 30
RMSE.C1=numeric(0)
RMSE.C5=numeric(0)

for (r in 1:30) {
set.seed(100+r)
l=sample(1:n,n*0.7)

dados <- data.frame(x,y)
dados.trein <- dados[ l,]
dados.teste <- dados[-l,]
x.teste <- dados.teste$x

# Ajustando modelo de regressão para C=1
mod.C1 <- ksvm(y~x,data=dados.trein,kernel='rbfdot',
            C=1,epislon=0.1)
y.est1 <- predict(mod.C1,dados.teste)
RMSE.C1[r]=sqrt(mean((dados.teste$y-y.est1)^2))

# Ajustando modelo de regressão para C=5
mod.C5 <- ksvm(y~x,data=dados.trein,kernel=rbfdot(),
               C=5,epislon=0.1)
y.est5 <- predict(mod.C5,dados.teste)
RMSE.C5[r]=sqrt(mean((dados.teste$y-y.est5)^2))
}

points(x.teste,y.est1,col='tomato',type='l',lwd=3)
points(x.teste,y.est5,col='#43a047',type='l',lwd=3)

legend(2,-2,c("C=1","C=5"),
       col=c('tomato','#43a047'),lty=1,
       lwd=3,box.col = "white")

#Sendo C=5 visualmente mais adequado que C=1
boxplot(RMSE.C1,RMSE.C5,ylab='RMSE',
        names=c("C=1","C=5"),col=c('tomato','#43a047'))

mean(RMSE.C5<RMSE.C1)

## [1] 0.8666667

Figura 6. Máquina de suporte de vetores para regressão

5.3 Outros casos

5.3.1 SVM para detecção de outliers

A ideia da utilização do SVM para detecção de outliers pode ser resumida através do cálculo de uma hiperesfera mínima que contenha todas as observações da base de dados considerados normais. A partir do parâmetro \(\nu\), é possível restringir ou relaxar o limite das margens, aumentando ou diminuindo a tolerância do que seria considerado um outlier. Mais detalhes sobre o modelo e suas particularidades podem ser encontrados em Schölkopf et al. (1999).

O exemplo de como utilizar o modelo pode ser visto abaixo.

# Gerando um novo conjunto de dados
set.seed(42)
n <- 1000
x <- rnorm(n = n,mean = 0,sd = 0.5)
n_out <- 5
sample_out <- sample(1:n,size = n_out)
x[sample_out] <- (rnorm(n = n_out,sd = 3))^2

# Gerando o modelo
library(kernlab)
outlier_svm <- ksvm(x,type= "one-svc",nu = 0.8, kernel = "vanilladot")

##  Setting default kernel parameters

# Visualizando os resultados
plot(1:n,x, type = "l")
outliers <- which(predict(outlier_svm))
points(outliers,x[outliers],pch=20, col = "#43a047",cex=3)

Figura 7. Máquina de suporte de vetores para detecção de outliers

5.3.2 SVM e Clustering

O SVM também pode ser utilizado para tarefas de agrupamento. Ben-Hur et al. (2001) desenvolveu o Support Vector Clustering (SVC), a partir da utilização do mapeamento das observações em espaços dimensionais de alta dimensão utilizando o kernel gaussiano. No espaço característico de alta dimensão, é então calculada a menor hiperesfera que contém \(\Phi(x_{i})\). Posteriormente, ao retornar ao espaço dimensional original, são obtidas as linhas de contorno que definem os limites de cada um dos clusters. Observações contidas em cada um desses limites pertencem ao mesmo cluster. À medida que o parâmetro de largura do kernel gaussiano é reduzido, o número de contornos desconexos no espaço de dados aumenta, levando a um número crescente de agrupamentos. Esse parâmetro é então responsável pela definição do número de clusters encontrados.

O SVC ainda pode lidar com outliers ao empregar os parâmetros de margem suave, permitindo que a esfera no espaço de características não inclua todos os pontos. Para valores grandes desse parâmetro, também é possível lidar com agrupamentos sobrepostos.

Figura 8. Exemplo de agrupamento de um conjunto de dados usando o SVC com \(C = 1\). Os vetores de suporte são designados por pequenos círculos, e as atribuições de cluster são representadas por diferentes cores pontos de dados. (a): \(q = 1\) (b): \(q = 20\) (c): \(q = 24\) (d): \(q = 48\). Fonte: Adaptado de Ben-Hur et al. (2001)

5.3.3 Análise de Sobrevivência

Os modelos de vetores de suporte para regressão também podem ser aplicados no contexto de análise de sobrevivência. Inicialmente, utilizar o SVR diretamente em uma análise de sobrevivência seria desafiador, visto que não é possível utilizar observações censuradas para estimar a função \(f(\mathbf{X})\). No entanto, Shivaswamy, Chu, e Jansche (2007) propuseram uma modificação ao tradicional SVR, de modo que pudesse adaptar a formulação original através de uma função de perda assimétrica que permite que observações censuradas (tanto à direita quanto à esquerda) pudessem ser processadas e computadas no problema de otimização inerente ao modelo. Dessa forma, o então originalmente nomeado Support Vector Regression for Censored Data (SVRc) - também conhecido como modelo de vetores de suporte para dados censurados - utiliza a alta capacidade de generalização do SVR, incluindo a alta flexibilidade em casos não lineares devido ao truque do kernel, em conjunto com os pressupostos necessários ao tipo de dado em questão.

A principal ideia está na modificação das Equações (4), (5) e (6) e as condições associadas ao problema de otimização. Inicialmente, define-se uma observação censurada como \((\mathbf{x}_{i},u_{i},l_{i})~ ~ ~\forall~i=1,\dots,n\), onde \(l_{i}\) é o limite inferior e \(u_{i}\) é o limite superior. No contexto de dados censurados, existem dois tipos diferentes de censura: os dados podem ser censurados à esquerda, com a tupla que representa a observação sendo dada por \((\mathbf{x}_{i},-\infty,u_{i})~ ~ ~\forall~i=1,\dots,n\), e os dados podem ser censurados à direita, com \((\mathbf{x}_{i},l_{i},+\infty) \forall~ ~i=1,\dots,n\). Se os dados não estão censurados, a observação é dada pela tupla \((\mathbf{x}_{i},y_{i}) \forall~~i=1,\dots,n\).

No SVRc, a função de perda, originalmente dada no SVR pela Equação (3), é agora definida por:

\[ L(f(\mathbf{x}_i),u_{i})=\max (0,l_{i}-f(\mathbf{x}_{i}),f(\mathbf{x}_{i})-u_{i})). \]

Seja \(I_u = \{i~|~l_i > -\infty , u_i < +\infty \}\), \(I_r = \{i~|~ u_i = +\infty \}\) e \(I_l = \{i~|~l_i = -\infty \}\). Então, \(I_u\) é o índice das amostras não censuradas, enquanto \(I_r\) e \(I_l\) são os índices das amostras que estão censuradas à direita e à esquerda, respectivamente. Definindo \(U=I_u \cup I_l\) e \(L=I_u \cup I_r\), a otimização da função objetivo convexa para o SVCR é então dada por:

\[ \min_{\mathbf{w},\mathcal{E},\mathcal{E}^{\ast}} \frac{1}{2} \mathbf{w} \boldsymbol{\cdot} \mathbf{w} + C\left(\sum_{ i \in U}\mathcal{E}_i+\sum_{i \in L}\mathcal{E}_i^{\ast}\right), \] \[ \text{sujeita à} \begin{cases} \hat{f}(\mathbf{x}_{i}) - u_{i} \leq \mathcal{E}_i^{\ast}, \hspace{0.2cm} \forall~ i \in U, \\ l_{i}-\hat{f}(\mathbf{x}_{i}) \leq \mathcal{E}_i, \hspace{0.2cm} \forall~ i \in L, \\ \mathcal{E}_i \geq 0 ~~\forall~ i \in U, \mathcal{E}_i^{\ast} \geq 0~~\forall~ i \in L. \end{cases} \] O restante de todas as expansões utilizando os multiplicadores Lagrangianos, a formulação dual, e a utilização do truque do kernel seguem naturalmente a partir das mesmas derivações das equações acima.

Um aspecto interessante do SVRc é acomodar diferentes representações em diferentes espaços dimensionais através do uso de diferentes funções de kernel em conjunto com o truque do kernel. Recentemente, Maia et al. (2023) mostraram como a exploração de diferentes funções de kernel, incluindo o wavelet kernel - não tão comum na literatura - pode melhorar de maneira significativa a capacidade preditiva do modelo em diferentes aplicações.

A implementação do SVRc está disponível em R através do pacote survivalsvm (Fouodo et al. (2018)). Um exemplo de sua utilização pode ser encontrado em detalhes na documentação do pacote, mas a demonstração abaixo segue utilizando o conjunto de dados veteran disponível no pacote survival (Therneau e Lumley (2015)).

# Para instalação dos pacotes:
# install.packages("survival")
# devtools::install_github("imbs-hl/survivalsvm")

# Carregando bibliotecas
library(survival)
library(survivalsvm)


# Criando os índices para amostra teste e treinamento
set.seed(42)
train_index <- sample(1:nrow(veteran),size = round(0.7*nrow(veteran)))
test_index <- (1:nrow(veteran))[-train_index]

# Criando o modelo
survivalsvm_reg <- survivalsvm(Surv(diagtime, status) ~ .,
                              subset = train_index, data = veteran,
                              type = "regression", gamma.mu = 1,
                              opt.meth = "quadprog", kernel = "add_kernel")

# Observando o resultado 
print(survivalsvm_reg)

## 
## survivalsvm result
## 
## Call:
## 
##  survivalsvm(Surv(diagtime, status) ~ ., subset = train_index, data = veteran, type = "regression", gamma.mu = 1, opt.meth = "quadprog", kernel = "add_kernel") 
## 
## Survival svm approach              : regression 
## Type of Kernel                     : add_kernel 
## Optimization solver used           : quadprog 
## Number of support vectors retained : 40 
## survivalsvm version                : 0.0.5

# Predizendo novas observações teríamos
pred.survsvm.reg <- predict(object = survivalsvm_reg,
                    newdata = veteran, subset = test_index)

# Acessando os cinco primeiros valores preditos
pred.survsvm.reg$predicted[1:5]

## [1] 11.69133 14.10848 11.49340 13.31303 13.03937

5.3.4 Séries Temporais

Devido a flexibilidade do SVR, muitos autores têm utilizado o SVR como um modelo preditivo para séries temporais (Pai et al. 2010; Makridakis, Spiliotis, e Assimakopoulos 2018). Esta situação não é adequada uma vez que o SVR é fundamentado pela teoria do aprendizado estatístico que supõe a independência das observações (V. Vapnik 1999). Assim, existe uma grande variedade de estudos recentes que incorporam a temporalidade na estrutura tradicional do SVR.

SVR Mediana é um método proposto por Hung, HUNG, e LIN (2014) o qual considerada uma suavização através de medianas para estrutura do SVR. Este método foi aplicado pelos autores em um estudo para previsão de saques em caixas eletrônicos na Inglaterra. A ideia gerral é que a mediana das médias móveis (MA - moving average) seja utilizada para previsão da variável \(Y_t,\) pela sua sequência de \(k\) médias aritméticas \(MA_k\). O uso da média aritmética pode ser visto como uma forma de suavização para séries temporais evitando flutuações indesejadas, porém, essa suavização pode também levar a média móvel a reagir aos eventos com algum atraso. Para exemplificar temos:

\[\begin{equation*} Y_t \approx\ MA_k = \frac{1}{k}(Y_{t-1} + Y_{t-2} + \dots + Y_{t-k+1}) \end{equation*}\]

No método foi aplicada a mediana de 4 semanas para modelar e prever a série temporal, a variável alvo utilizada é a diferença entre o valor real e a mediana semanal da média móvel. Essa diferença vai ser chamada de \(D_t\) que vai ser a variável resposta a ser predita pelo SVR.

\[\begin{equation*} \begin{gathered} \text{week}_{MA_2} = \frac{1}{2}\biggl(Y_{t-7} + Y_{t-14}\biggr) \\ \text{week}_{MA_3} = \frac{1}{3}\biggl(Y_{t-7} + Y_{t-14} + Y_{t-21}\biggr) \\ \text{week}_{MA_4} = \frac{1}{4}\biggl(Y_{t-7} + Y_{t-14} + Y_{t-21} + Y_{t-28}\biggr) \\ \\ Y_t \approx D_t +\text{mediana}\biggl(\text{week}_{MA_2}, \text{week}_{MA_3}, \text{week}_{MA_4}\biggr) \end{gathered} \end{equation*}\]

Além disso, existem outros métodos para a incorporação da temporalidade no SVR. O SRV-Recorrente (Chen, Jeong, e Härdle 2008) utiliza uma estrutura de recorrência baseada nos modelos ARMA, bem como o método SRV-Recorrente com estrutura GARCH (Chen, Härdle, e Jeong 2010) que trata da recorrência, porém se baseando na modelagem da volatilidade. Existem também propostas mais recentes de modelos híbridos, como o modelo SVR-wavelet (Raimundo e Okamoto 2018), um modelo de predição adaptativo e híbrido que integra modelos wavelet e Support Vector Regression (SVR), para predição de séries temporais.

Ainda não há implementações disponíveis em pacotes para SVR-Mediana, SVR-Recorrente e SVR-GARCH. Porém, SVR-wavelet para séries temporais está disponível em Linguagem R através do pacote WaveletSVR. O script abaixo exemplifica o uso SVR-wavelet para uma série temporal gerada através de um modelo ARMA(2,2).

require('WaveletSVR')

#Geração dos dados
set.seed(100)

data=arima.sim(n = 100, list(ar = c(0.9, -0.5), ma = c(-0.20, 0.25)),
         sd = 1)

#Ajuste do modelo
WSVR<-WaveletFittingsvr(ts=data,tlag=2,Waveletlevels=3)

## Gráfico da série temporal
plot(data,type='l',lwd=2,col='gray')
lines(3:80,WSVR$TrainFittedValue,col='red',lty=2,lwd=1.5)
lines(81:100,WSVR$TestPredictedValue,col='#2e7d32',lty=2,lwd=1.5)
legend(0,-2.3, c("treino","teste"),col=c("red","#2e7d32"),lwd=2,lty=3,cex=0.5)