SAUDAÇÕES!

Seja bem vindo à página do professor Pedro Albuquerque. Para saber mais sobre meu currículo, disciplinas ministradas e interesses de pesquisa, navegue no menu disponível no topo da página.

sexta-feira, 15 de setembro de 2017

CAViaR - Conditional Autoregressive Value at Risk by Regression Quantiles

O modelo CAViaR - Conditional Autoregressive Value at Risk by Regression Quantiles é um modelo utilizado para se mensurar o Value at Risk para períodos diários e, assim, monitorar o risco de maneira autogregressiva.

Engle e Manganelli (2004) sugerem quatro possíveis funções para a estimação do CAViaR:

  • Adaptive: $VaR_{t}=VaR_{t-1}+\beta_{1}\left\{\left[1+\exp\left(G\left[y_{t-1}-VaR_{t-1}\right]\right)^{-1}\right]-\tau\right\}$.
  • Symmetric absolute value: $VaR_{t}=\beta_{1}+\beta_{2}VaR_{t-1}+\beta_{3}|y_{t-1}|$.
  • Asymmetric slope: $VaR_{t}=\beta_{1}+\beta_{2}VaR_{t-1}+\beta_{3}(y_{t-1})^{+}+\beta_{4}(y_{t-1})^{-}$.
  • Indirect GARCH(1, 1): $VaR_{t}=\left(\beta_{1}+\beta_{2}VaR_{t-1}^{2}+\beta_{3}y_{t-1}^{2}\right)^{1/2}$

onde $y_{t}$ é o retorno do ativo no tempo $t$ com $t=1,\dots, T$. $G$ é uma constante positiva (Engle e Manganelli (2004) sugerem $G=10$), $\tau$ é o nível do VaR desejado (usualmente $\tau=0.01$ ou $\tau=0.05$) e $(x)^{+}=\max(x,0)$ e $(x)^{-}=-\min(x,0)$

Como as funções são recursvivas é necessário uma primeira estimativa para o $VaR_{t=1}$. Considerando os mesmos dados de Engle e Manganelli (2004) os quais trabalharam com os retornos dos preços das açoes da GM, IBM e SP500 para Abril 7, 1986 até Abril 7, 1999, podemos obter a primeira estimativa de $VaR_{t=1}$ como:

#Limpa os workspace
rm(list=ls())
#Lê os dados
stocks<-read.table("dataCAViaR.txt")
colnames(stocks)<-c("GM","IBM","SP500")
#Nível do VaR desejado
tau<-0.05
#Divide a base em treinamento e validação
train<-stocks[1:2892,]
valid<-stocks[2893:3392,]
#Criar o vetor de VaR
VaR<-rep(NA,nrow(train))
#Primeira estimativa
VaR[1]<- -qnorm(tau)*sd(train$GM)
Em seguida, Engle e Manganelli (2004) sugerem estimar os parâmetros dos modelos por meio da estrutura de Regressão Quantílica: $\min_{\boldsymbol\beta}\frac{1}{T}\sum_{t=1}^{T}\left[\tau-1\left(y_{t} \leq VaR_{t}\right)\right]\left[y_{t}-VaR_{t}\right]$ onde $1(\cdot)$ é uma função indicadora que assume valor igual a 1 quando seu argumento é verdadeiro $y_{t}\lt VaR_{t}$ e 0 caso contrário. Assim, a função objetivo para o modelo Adaptive é dada por:
#Sugestão Engle e Manganelli (2004)
G<-10
#Function
adpatative<-function(x){
  beta1<-x[1]
  for(i in 2:nrow(train)){
    VaR[i]<- VaR[i-1] + beta1*((1+(exp(G*(train$GM[i-1]-VaR[i-1]))^(-1)))-tau)
  }
  #Objective Function
  res<-sum((tau-(train$GM < VaR))*(train$GM-VaR))/nrow(train)
  return(res)
}
O próximo passo é otimizar essa função, isso pode ser feito utilizando-se a função optim:
#Realiza a otimização usando o métood de Brent
res<-optim(par=c(0), adpatative, method="Brent", lower=c(-10), upper=c(10))
#Resultado do parâmetro de interesse
beta<-res$par
A função optim tem como argumentos nesse exemplo: valor inicial para o(s) parâmetro(s) de interesse, função a ser minimizada, método utilizado e limites inferiores e superiores para o(s) parâmetro(s) de interesse. O último passo é a previsão do VaR:
VaR<-rep(NA,nrow(dados))
VaR[1]<- -qnorm(tau)*sd(dados$GM)
#Previsão do CAViaR
adpatativeForecast<-function(x){
  beta1<-x[1]
  for(i in 2:nrow(dados)){
    VaR[i]<- VaR[i-1] + beta1*((1+(exp(G*(dados$GM[i-1]-VaR[i-1]))^(-1)))-tau)
  }
  #Objective Function
  return(VaR)
}
onde o objeto dados pode ser as bases de treinamento e validação:
#Base de treinamento
dados<-train
#Faz a previsão
VaR.train<-adpatativeForecast(beta)
#Medida de acuracia hits
mean(dados$GM< -VaR.train)-tau

O mesmo processo pode ser realizado para as demais funções.

terça-feira, 15 de agosto de 2017

RcppEigen: Operações Matriciais e Rcpp - Parte 3

A biblioteca RcppEigen é uma das mais eficientes para o cômputo numérico de operações algébricas uma vez que utiliza algortimos optimizados para a memória cache e SIMD.

Por exemplo, considere o desafio de estimar uma Matriz de Variâncias e Covariâncias com base em um conjunto de dados. Essa estimação levaria em conta pelo menos $N\times(N-1)/2$ operações (calculando todas as variâncias e covariâncias na matriz triangular uma vez que a matriz é simétrica). Esse procedimento levaria em consideração o uso de loops aninhados (nested loops) os quais podem ser muito ineficientes e demorados quando implementados diretamente em R.

Nesse sentido, podemos usar as operações do pacote RcppEigen também documentadas na página do Eigen C++:

// [[Rcpp::export]]
#include 
// [[Rcpp::depends(RcppEigen)]]
#include 
using namespace Rcpp;

// [[Rcpp::export]]
Eigen::MatrixXd covCalcula(Eigen::MatrixXd mat) {
  //Centra as variáveis com respeito a média
  Eigen::MatrixXd centered = mat.rowwise() - mat.colwise().mean();
  //Calcula a covariância
  Eigen::MatrixXd cov = (centered.adjoint() * centered) / double(mat.rows() - 1);
  return cov;
}

Observe que as operações .rowwise() e .colwise() se referem, respectivamente, aos cômputos para cada linha e para cada coluna. Já o método .adjoint() calcula a Matriz Adjunta.

A invocação desse método no R é realizada da seguinte forma:

#Limpa o Working Directory
#Chama um conjunto de dados
data("ChickWeight")
#Converte todas as variáveis para numeric
dat2 <- data.frame(lapply(ChickWeight, function(x) as.numeric(as.character(x))))
#Transforma para matriz
mat<-as.matrix(dat2)
#Calcula a matriz de covariância
covCalcula(mat)


segunda-feira, 17 de julho de 2017

RcppEigen: Operações Matriciais e Rcpp - Parte 2


Uma das vantagens de se trabalhar com o RcppEigen e Rcpp, além da velocidade, é o fato de que não é necessário implementar todas as funções de interesse em C++. Podemos importar funções já existentes do R para que essas sejam executadas no ambiente Rcpp.

Por exemplo, um problema comum na formação de portfólios é o fato das matrizes de correlações estimadas entre os retornos dos ativos não ser Positivas Definidas. Nesses casos, problemas teóricos e numéricos podem surgir impedindo a obtenção de uma solução ótima para o problema de alocação de portfólios. Higham, N. (2002) propôs uma abordagem numérica para se encontrar a matriz Positiva Definida mais próxima da matriz de ineteresse.

Essa função está disponível no pacote Matrix e é denominada nearPD. Podemos importá-la no ambiente Rcpp da seguinte forma:

// [[Rcpp::export]]
Eigen::MatrixXd nearPDefinite(Eigen::MatrixXd X, int maxit=1e+6, double eigtol = 1e-06, double conv_tol = 1e-07, double posd_tol = 1e-08, bool keepDiagonal=false){
  //Cria o objeto com o Environment.
  Rcpp::Environment G = Rcpp::Environment::global_env();
  //Objeto com as informações do pacote Matrix
  Rcpp::Environment Matrix("package:Matrix");
  //Importa a função de ineteresse
  Function f = Matrix["nearPD"];
  //Obtêm o resultado da função nearPD
  Rcpp::List res = f(X, false, keepDiagonal, true, false, true, false, false,  eigtol, conv_tol, posd_tol, 100, "I", false);
  //Pega o primeiro elemento da lista e transforma para matrix
  Rcpp::NumericMatrix mat = internal::convert_using_rfunction(wrap(res[0]), "as.matrix");
  //Converte de Rcpp::NumericMatrix para Eigen::MatrixXd
  Eigen::MatrixXd Xpd=convertMatrix(mat);
  //Retorna o resultado
  return(Xpd);
}

A função anterior, após compilada, pode ser utilizada no ambiente R de maneira direta pelo nome nearPDefinite. Entretanto, é necessário que a biblioteca Matrix tenha sido habilitada e esteja disponível no environment do R para que seja utilizada. Essa estratégia pode ser para outras funções de outros pacotes também.

quinta-feira, 15 de junho de 2017

RcppEigen: Operações Matriciais e Rcpp - Parte 1


Uma das grandes críticas que os usuários do R fazem é que por vezes algumas operações são muito demoradas no R. Em especial, a execução de loops.

Há, no entanto alguns truques que podem ser utilizados para se aumentar a velocidade dos códigos produzidos em R, quais sejam: usar operações algébricas vetoriais ao invés de loops, utilizar os membros da família apply, utilizar processamento em paralelo , etc.

Entretanto, algumas vezes é necessário programar funções mais complexas e que necessitam do uso de loops clássicos. Quando esse for o caso, uma boa solução é utilizar o pacote RcppEigen que é uma versão em R da biblioteca em C++ Eigen.

Nesse grupo de posts acerca do RcppEigen vamos mostrar como esse pacote pode ser usado e também como novos pacotes para o R podem ser criados. O primeiro passo é construir um projeto de pacote no RStudio:


Definido a criação de um novo pacote o próximo passo é determinar:


Podemos criar um projeto com base em um diretório já existente ou em um novo diretório. Escolhendo um novo diretório temos as seguintes opções:


Escolhemos nesse caso a opção R Package. Em seguida, definimos o local e nome do pacote que desejamos criar:


Automaticamente o R cria alguns arquivos e pastas. Para usar o RcppEigen precisamos de alguns arquivos específicos. Esses arquivos são obtidos executando uma função padrão do pacote:

#Limpa o Workspace
#Habilita o pacote
library(RcppEigen)
#Define o local temporário para criação dos arquivos
setwd("C:\\LocalTemorario")
#Cria os arquivos necessários
RcppEigen.package.skeleton("ArquivosRcppEigen")

Na pasta definida pelo comando setwd um conjunto de arquivos é criado.

Deve-se copiar os arquivos DESCRIPTION e NAMESPACE e também a pasta src (todas elas dentro da pasta "ArquivosRcppEigen") para a pasta do pacote criado no RStudio (nesse exemplo denominada "NomeDoPacote"). Após essas etapas o pacote estará pronto para receber as funções escritas em RcppEigen, as quais veremos nas próximas partes desse tutorial.

segunda-feira, 15 de maio de 2017

taskscheduleR: Agendando tarefas no R.


A automação de tarefas é importante quando deseja-se executar um mesmo script de maneira regular. Nesse sentido o pacote taskscheduleR pode apoiar no agendamento de scripts feitos em R.

Como exemplo vamos supor que estamos interessados em decidir se devemos comprar ou vender uma ação no dia seguinte baseado no preço de fechamento dos dias anteriores usando um modelo de Máquinas de Suporte Vetorial.

#Limpa o Workspace
rm(list=ls())
#Chama a biblioteca
library(devtools)
library(quantmod)
library(kernlab)
#Instala o pacote taskscheduleR
devtools::install_github("jwijffels/taskscheduleR")
#Chama a biblioteca taskscheduleR
library(taskscheduleR)
#Especifica as datas de interesse
startDate = as.Date("2013-01-01") 
endDate = as.Date("2015-12-31")
#Obtêm os dados do ativo PETR4 e PETR3
getSymbols("PETR4.SA", src="yahoo",from=startDate,to=endDate)
#Faz o Candle Chart
candleChart(PETR4.SA,theme='white', type='candles') 
reChart(major.ticks='months',subset='first 16 weeks') 

O Candle Chart é apresentado a seguir:


Em seguida vamos treinar uma máquina para tentar apontar a direção do resultado futuro, isto é, se a ação irá apresentar um resultado positivo ou negativo:

#Retorno do Preço de Fechamento
ret<-na.omit(ClCl(PETR4.SA))
#Direção do mercado
y<-as.numeric(ifelse(ret<0,0,1))
#Cria a base de dados
df<-data.frame("Y"=c(y[-1],NA),"X"=ret)
#Renomeia as colunas
colnames(df)<-c("Y","X")
Nesse bloco a ideia é tentar prever a direção $Y=sinal(r_{t})$ usando o retorno do período anterior $X=r_{t-1}$. Supondo que o processo de treinamento tenha sido realizado e que tenhamos descobertos $C=0.1$ e $\sigma=16$ para o modelo C-SVC e kernel Gaussiano, a máquina é construída fazendo-se:
#Constrói a máquina
svm <- ksvm(Y~X,data=df, type="C-svc",C=0.1,
            kpar=list(sigma=16),cross=3)
#Faz a previsão
pred<- kernlab::predict(svm,df)
#Compara com os valores reais
table(pred,df$Y)
#Salva a máquina treinada
save.image("Maquina.RData")
Supondo que essa acurácia seja a desejada, agora vamos programar para que o R todo dia baixe a cotação do dia anterior e faça a previsão de qual deveria ser nossa ação no dia atual (vender se a previsão é que o retorno seja positivo ou comprar se a previsão seja do retorno negativo). As opções do pacote taskscheduleR são: ‘ONCE’, ‘MONTHLY’, ‘WEEKLY’, ‘DAILY’, ‘HOURLY’, ‘MINUTE’, ‘ONLOGON’, ‘ONIDLE’. O primeiro passo é montar um script geral que será executado todo dia:
#Chama a biblioteca
library(quantmod)
library(kernlab)
#Lê a máquina criada
load("Maquina.RData")
#Obtêm a cotação atual
getSymbols("PETR4.SA", src="yahoo")
#Último retorno observado
ret<-tail(ClCl(PETR4.SA),1)
#Faz a previsao para o dia
previsao<-kernlab::predict(svm,ret)
#Monta o arquivo de saida
fileConn<-file("Decisao.txt")
writeLines(ifelse(previsao==1,"Vender","Comprar"), fileConn)
close(fileConn)
Esse script deve ser salvo em um arquivo de extensão .R (por exemplo, ScriptPrevisao.R) e, no nosso exemplo, ele deverá ser executado todo dia as 8:00 AM. Após a execução, ele criará um arquivo denominado Decisao.txt que contêm a decisão do dia: Comprar ou Vender. O último passo é programar o computador para executar esse script todo dia no horário determinado, para isso, basta executarmos uma única vez no R:
#Faz o agendamento
library(taskscheduleR)
myscript <- system.file("saida", "ScriptPrevisao.R", package = "taskscheduleR")
#O computador irá executar todo dia as 8:00AM
taskscheduler_create(taskname = "DecisaoPETR4", rscript = myscript, 
                     schedule = "DAILY", starttime = "08:00")
Caso queira deletar o agendamento para que ele não se repita mais, basta fazer:
#Deleta o agendamento
library(taskscheduleR)
taskscheduler_delete(taskname = "DecisaoPETR4")




sábado, 15 de abril de 2017

Introdução ao RCpp - Parte 4.


Frequentemente é necessário inicializar um NumericVector com algum valor de constante fixo, nesse caso podemos usar a função fill(.):

#include 
using namespace Rcpp;
// Função para repetir o valor no vetor
// [[Rcpp::export]]
NumericVector inicializaVetor(double valor, int tamanho) {
  NumericVector x(tamanho);
  x.fill(valor);
  return x;
}

Após compilar o código, a invocação no R é quase imediata:

#Invoca a biblioteca RCpp
library(Rcpp)
#Define o endereço do arquivo Cpp
setwd("C:\\Blog\\Source")
#Compila o código em C++
sourceCpp('Teste.cpp')
inicializaVetor(2,5)

O pacote Rcpp tem muitas funcionalidades. Por meio do Rcpp podemos por exemplo, usar funções já implementadas no R, como:

// [[Rcpp::export]]
NumericVector calculaCurtose(NumericVector x) {
  Rcpp::Environment funcoesPacote("package:e1071");
  Rcpp::Function funcaoCurtose = funcoesPacote["kurtosis"];
  return funcaoCurtose (x);
}

Ou por exemplo, gerar números com a distribuição normal:

// [[Rcpp::export]]
NumericVector geraNormal(){
  Rcpp::Environment stats("package:stats");
  Rcpp::Function rnorm = stats["rnorm"];
  return rnorm(10, Rcpp::Named("sd", 100.0));
}

Comandos esses que após compilados são invocados na forma:

#Invoca a biblioteca RCpp
library(Rcpp)
#Define o endereço do arquivo Cpp
setwd("C:\\Blog\\Source")
#Compila o código em C++
sourceCpp('Teste.cpp')
#Executa as funções:
x<-rnorm(1919)
calculaCurtose(x)
geraNormal()
Entre os objetos mais flexíveis do Rcpp está o objeto GenericVector o qual funciona como uma lista permitindo a inserção de qualquer outro objeto em seus slots. Além do mais, é possível usar os seguintes métodos nesse objeto:
  • push_back(x) - Insere o elemento x no final da lista.
  • push_front(x) - Insere o elemento x no início da lista.
  • insert(i, x) - Insere o elemento x na posição i.
  • erase(i) - Apaga o elemento existente na posição i.

Todas essas funções são dinâmicas e não há a necessidade de criar uma nova lista para que os novos elementos sejam inseridos.

quarta-feira, 15 de março de 2017

Introdução ao RCpp - Parte 3.


Aqui neste post continuaremos com a análise dos possíveis métodos desenvovidos usando Rcpp. O primeiro exemplo trata do uso de matrizes, especificamente, vamos criar uma função para fazer a transposição de uma matriz:

#include 
using namespace Rcpp;
// Função para transposição de matriz
// [[Rcpp::export]]
NumericMatrix transpose(NumericMatrix x) {
  //Obtêm o número de linhas e colunas:
  int nrow = x.nrow(), ncol = x.ncol();
  //Cria a matriz com resultado
  NumericMatrix resultado(ncol,nrow);
  //Faz a transposição:
  for(int i=0;i < nrow; i++){
    for(int j=0;j < ncol; j++){
      resultado(j,i)=x(i,j);
    }
  }
  return resultado;
}
Após compilar o código, a invocação no R é quase imediata:
#Invoca a biblioteca RCpp
library(Rcpp)
#Define o endereço do arquivo Cpp
setwd("C:\\Blog\\Source")
#Compila o código em C++
sourceCpp('Teste.cpp')
#Cria uma matriz de exemplo
mat<-matrix(c(2,8,9,0,0,0,4,3,2), nrow=3,ncol=3)
#Transpoem a matriz
transpose(mat)
Métodos mais sofisticados usando funções implementadas em R e invocadas em C++ também podem ser utilizadas:
//Usando funções do R:
// [[Rcpp::export]]
NumericVector ExemploKernel(NumericVector x, NumericVector y, Function f) {
  NumericVector res = f(x,y);
  return res;
}
Método esse que recebe como argumento dois vetores numéricos e uma função bivariada:
#Cria uma função no R
gaussian<-function(x,y){
  return(sum((x-y)^2))
}
#Gera dois vetores
x<-rnorm(1000)
y<-rnorm(1000)
#Invoca a função em Rcpp:
ExemploKernel(x,y,gaussian)