Wprowadzenie

Strategia pair trading działa poprzez wykorzystywanie krótkookresowych odchyleń między cenami akcji, które w długim okresie zachowują równowagę. W momencie, gdy taka krótkookresowa nierównowaga ma miejsce, jedna ze spółek staje się relatywnie przeszacowana w stosunku do drugiej. Można wtedy zainwestować w portfel składający się z takiej pary, w którym zajmujemy krótką pozycję w spółce relatywnie przeszacowanej i długą w relatywnie niedoszacowanej. Pozycję zamykamy poprzez zajęcie przeciwnych pozycji w momencie, kiedy ceny spółek wrócą do poziomu równowagi długookresowej. Zgodnie z założeniami analizie poddano pięć par skointegrowanych spółek, w okresach in-sample (01-01-2011 - 31-12-2013) oraz out-of-sample (01-01-2014 - 30-06-2015). W celu zbadania stopnia zintegrowania zmiennych użyto testu Dickeya-Fullera. Opracowano sygnały otwarcia/zamknięcia pozycji na podstawie wartości nierównowagi krótkookresowej wyznaczonej przez reszty z modelu. Opracowano również wykresy equity line pokazujące zyskowność analizowanej strategii w porównaniu z portfelem opartym na prostej stratgii buy&hold dla indeksu WIG20.

Pierwsza para spółek - Enea i PGNiG

Okres in-sample

Wykorzystałem dzienne notowania dla spółek Enea oraz PGNiG. Dane pochodzą z serwisu http://stooq.com.

1. DANE

Ładuję do pamięci wymagane pakiety. Importuję dane dotyczące pierwszej spółki oraz dokonuję niezbędnych poprawek: poprawiam format daty oraz usuwam niepotrzebne zmienne:

library(xts)
library(PerformanceAnalytics)
library(tseries)
url1 <- "http://stooq.com/q/d/l/?s=ena&i=d"
pkn <- read.csv(url1,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
pkn$Date <- as.Date(pkn$Date)
pkn <- pkn[, c("Date","Open", "Close")]

Podobnie postępuję z drugą spółką oraz notowaniami WIG20:

url2 <- "http://stooq.com/q/d/l/?s=pgn&i=d"
lts <- read.csv(url2,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
lts$Date <- as.Date(lts$Date)
lts <- lts[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20 <- read.csv(url1,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
wig20$Date <- as.Date(wig20$Date)
wig20 <- wig20[, c("Date","Open", "Close")]

Ograniczam zakres notowań do okresu in-sample oraz przekształcam datę na właściwy format:

pkn<-pkn[!(pkn$Date < "2011-01-01"), ]
pkn<-pkn[!(pkn$Date > "2013-12-31"), ]

lts<-lts[!(lts$Date < "2011-01-01"), ]
lts<-lts[!(lts$Date > "2013-12-31"), ]

wig20<-wig20[!(wig20$Date < "2011-01-01"), ]
wig20<-wig20[!(wig20$Date > "2013-12-31"), ]

pkn$Date <- as.Date(pkn$Date, format = "%Y-%m-%d")
lts$Date <- as.Date(lts$Date, format = "%Y-%m-%d")
wig20$Date <- as.Date(wig20$Date, format = "%Y-%m-%d")

Wykres cen zamknięcia:

plot(pkn$Date, pkn$Close,type="l", ylim=c(0, 25), ylab = "Cena zamknięcia", xlab = "Data", col="lightblue3")
lines(pkn$Date, lts$Close, type="l", col="firebrick3")
title("Ceny zamknięcia Enea i PGNiG (in-sample)", cex.main=0.85)

2. Obliczam pierwsze różnice:

lts$dlts <- diff.xts(lts$Close)
pkn$dpkn <- diff.xts(pkn$Close)

W analizie wykorzystam funkcję testdf z zajęć.

testdf <- function(variable, adf_order){
  library(urca)
  library(lmtest)
  results_adf <- data.frame(order = -1, adf = 0,
                            p_adf = "", bgodfrey = 0, p_bg = 0)
  variable <- variable[!is.na(variable)]
  
  for(order in 0:adf_order){
    df.test_ <- ur.df(variable, type = c("drift"), lags = order)
    df_ <- df.test_@teststat[1]
    df_crit <- df.test_@cval[1, ]
    df_crit <- (df_ < df_crit) * 1
    p_adf <- ifelse (sum(df_crit) == 0,
                     ">10pct",
                     paste("<", names(df_crit)[min(which(df_crit == 1))],
                           sep = ""))
    resids_ <- df.test_@testreg$residuals
    bgtest_ <- bgtest(resids_ ~ 1, order = 1)
    bgodfrey <- bgtest_$statistic
    names(bgodfrey)<-NULL
    p_bg <- bgtest_$p.value
    
    results_adf <- rbind(results_adf,
                         data.frame(order = order,
                                    adf = df_,
                                    p_adf = p_adf,
                                    bgodfrey = bgodfrey,
                                    p_bg = p_bg)
    )
  }
  
  results_adf<-results_adf[results_adf$order>=0,]
  
  plot(variable,
       type = "l",
       col = "blue",
       lwd = 2,
       main = "Plot of the examined variable")
  
  return(results_adf)
}

3. Sprawdzam stopień zintegrowania zmiennych

testdf(variable = pkn$Close,  adf_order = 1)

##   order       adf p_adf    bgodfrey      p_bg
## 2     0 -2.958826 <5pct 1.087479435 0.2970304
## 3     1 -2.918316 <5pct 0.004119579 0.9488238
testdf(variable = pkn$dpkn, adf_order = 1)

##   order       adf p_adf    bgodfrey      p_bg
## 2     0 -28.44080 <1pct 0.006012542 0.9381935
## 3     1 -20.98082 <1pct 0.007458591 0.9311777
testdf(variable = lts$Close,  adf_order = 1)

##   order       adf  p_adf    bgodfrey      p_bg
## 2     0 -1.382862 >10pct 3.751289060 0.0527668
## 3     1 -1.314990 >10pct 0.003868325 0.9504069
testdf(variable = lts$dlts, adf_order = 1)

##   order       adf p_adf    bgodfrey      p_bg
## 2     0 -29.32265 <1pct 0.004374079 0.9472689
## 3     1 -20.59441 <1pct 0.001916712 0.9650796

Statystyka testu Dickeya-Fullera dla spółki Enea wynosi -2.958826. Oznacza to brak podstaw do odrzucenia hipotezy zerowej o braku stacjonarności zmiennej losowej. Następnie przeprowadziłem test dla pierwszych róznic cen akcji spółki Enea. Wartość statystykli DF wyniosła -28.44080. Oznacza to odrzucenie hipotezy zerowej, że zmienna jest zintegrowana co najmniej w stopniu 2. Wartość statystyki Breuscha-Godfreya równa 0.006012542 oraz p-value wynoszące 0.9381935 świadczą o tym, że realizacja jednokrotnie zróżnicowanej zmiennej jest białym szumem. Świadczy o tym również kształt otrzymanego wykresu.

Wynika z tego, że ceny akcji Enea są zintegrowane ~I(1).

Odpowiednie wartości statystyk dla akcji spółki PGNiG wyniosły:

  1. -1.382862
  2. -29.32265
  3. 0.004374079
  4. 0.9472689

Na podstawie powyższych wartości jak również analizy otrzymanego wykresu stwierdzam, że ceny akcji PGNiG są zintegrowane ~I(1)

Obie zmienne są ~I(1) więc sprawdzam, czy są skointegrowane.

4. Szacuję wektor kointegrujący

model.coint <- lm(pkn$Close ~ lts$Close)
summary(model.coint)
## 
## Call:
## lm(formula = pkn$Close ~ lts$Close)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -3.3067 -1.0166  0.0070  0.7617  4.0743 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 20.72018    0.25242   82.08   <2e-16 ***
## lts$Close   -1.37486    0.05824  -23.61   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.383 on 745 degrees of freedom
## Multiple R-squared:  0.4279, Adjusted R-squared:  0.4272 
## F-statistic: 557.3 on 1 and 745 DF,  p-value: < 2.2e-16

5. Testuję stacjonarność reszt

testdf(variable = model.coint$residuals, adf_order = 1)

##   order       adf p_adf     bgodfrey      p_bg
## 2     0 -3.720892 <1pct 0.6407312708 0.4234461
## 3     1 -3.639094 <1pct 0.0001807894 0.9892721

W celu zbadania kointegracji obu zmiennych należy oszacować wektor kointegrujący i przetestować stacjonarność reszt.

Wartość statystyki DF równa -3.720892 oznacza brak podstaw do odrzucenia hipotezy zerowej o stacjonarnosci reszt.

Wartość statystyki BG równa 0.6407312708 oraz wartość p-value na poziomie 0.4234461 oznaczają brak podstaw do odrzucenia hipotezy o braku autokorelacji reszt.

Oznacza to, że reszty z oszacowanego modelu są stacjonarne, co świadczy o kointegracji spółek Enea i PGNiG.

6. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Skoro ceny obu spółek są skointegrowane, będę mógł zastosować strategię pair-trading. W tym celu muszę wyznaczyć sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd

model.coint$czy_sd<-0
model.coint$czy_zero<-0
for(i in 1:747){
  model.coint$czy_sd[i] = ifelse(abs(model.coint$residuals[i])>sd(model.coint$residuals),1,0)
}

Sprawdzam, czy dwie kolejne wartości przecięły wartość 0.

for(i in 1:746){
  model.coint$czy_zero[i+1]=ifelse(model.coint$residuals[i]*model.coint$residuals[i+1]>0,0,1)
}

model.coint$zamk<-0
model.coint$otw<-0
i=1
k=1

for(k in 1:length(model.coint$residuals)){
  model.coint$zamk[k]<-NA
}
for(k in 1:length(model.coint$residuals)){
  model.coint$otw[k]<-NA
}

Uzupełniam kolumny otw oraz zamk sygnałami.

while(i <= 747){
  if(model.coint$czy_sd[i]==1){
    model.coint$otw[i]=ifelse(model.coint$residuals[i]>sd(model.coint$residuals),sd(model.coint$residuals),-sd(model.coint$residuals))
    i<-i+1
    while(i <= 747){
      if(model.coint$czy_zero[i]==1){
        model.coint$zamk[i]=0
        break
        }
      i<-i+1
    }
  }
  i<-i+1
}

7. Wykres pozycji dla obu spółek

plot(pkn$Date, model.coint$residuals, type="l", ylab = "Nierównowaga", xlab = "Data")
abline(h=c(sd(model.coint$residuals),-sd(model.coint$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pkn$Date, y=model.coint$otw, col="red", pch=17)
points(x=pkn$Date, y=model.coint$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (in-sample)", cex.main=0.85)

Powyższy wykres pokazuje sygnały otwarcia/zamknięcia pozycji w strategii pair-traiding w okresie in-sample. Opierają się one na wykresie wartości nierównowagi krótkookresowej, czyli tak naprawdę na wartościach reszt z oszacowanego modelu. Moment wychylenia się nierównowagi krótkookresowej poza wartość jednego odchylenia standardowego (czerwona linia) stanowi sygnał otwarcia pozycji (czerwony trójkąt). Pozycję zamykam (niebieski kwadrat), kiedy wartość reszty wróci do średniej, czyli wartości 0 (niebieska linia).

8. Obliczam zwroty dla obu spółek.

Zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

i=1
for(i in 2:747){
  model.coint$zwrot_A[i] = pkn$Open[i]/pkn$Open[i-1]
}

i=1
for(i in 2:747){
  model.coint$zwrot_B[i] = lts$Open[i]/lts$Open[i-1]
}

model.coint$otw[model.coint$otw==sd(model.coint$residuals)]<-1
model.coint$otw[model.coint$otw==-sd(model.coint$residuals)]<-1
model.coint$otw[is.na(model.coint$otw)]<-0
model.coint$zamk[model.coint$zamk==0]<-1
model.coint$zamk[is.na(model.coint$zamk)]<-0

9. Tworzę kolumnę pozycja.

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1
while(i<=746){
  repeat{
    model.coint$pozycja[i]=0
    i<-i+1
    if(model.coint$otw[i-1]!=0||i>746){break}
  }
  while(model.coint$zamk[i]!=1&&i<=746){
    model.coint$pozycja[i]=ifelse(model.coint$residuals[i]>0,1,-1)
    i<-i+1
  }
  model.coint$pozycja[i]=ifelse(model.coint$residuals[i-1]>0,1,-1)
  i<-i+1
}

10. Obliczam zwrot z pozycji.

Kolumna zwrot z pozycji mówi, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
model.coint$zwrot_z_pozycji[1]=0
for(i in 2:747){
  model.coint$zwrot_z_pozycji[i] = model.coint$pozycja[i]*(model.coint$zwrot_B[i]-model.coint$zwrot_A[i])
}

i=2
for(i in 2:747){
  if(model.coint$otw[i]==1){
    model.coint$zwrot_z_pozycji[i+1] = model.coint$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:747){
  if(model.coint$zamk[i]==1){
    model.coint$zwrot_z_pozycji[i] = model.coint$zwrot_z_pozycji[i]*0.999
  }
}

11. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

model.coint$portfel[1]=1
i=2
for(i in 2:747){
  model.coint$portfel[i] = model.coint$zwrot_z_pozycji[i-1]*model.coint$portfel[i-1]+model.coint$portfel[i-1]
}

12. Obliczam dzienne zwroty z indeksu WIG20

Liczę też zwrot z portfela dla strategii buy&hold dla WIG20.

i=2
wig20$zwrot<-0
for(i in 2:747){
  wig20$zwrot[i] = wig20$Open[i]/wig20$Open[i-1]
}

wig20$portfel[1]=1
i=2
for(i in 2:747){
  wig20$portfel[i] = wig20$zwrot[i]*wig20$portfel[i-1]
}

13. Wykres equity line

plot(pkn$Date, model.coint$portfel, type="l", ylim=c(0.5, 2), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pkn$Date, wig20$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "topleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (in-sample)", cex.main=.85)

Praktycznie do połowy roku 2013 widoczna jest wyraźna przewaga zyskowności strategii opartej na kointegracji nad buy&hold na WIG20. Z czasem jednak jej zyskowność słabła spadając pod koniec analizowanego okresu poniżej poziomu rentowności strategii buy&hold.

Okres out-of-sample

1. DANE

Importuję dane i dokonuję niezbędnych poprawek na tych samych spółkach w okresie out-of-sample.

url1 <- "http://stooq.com/q/d/l/?s=ena&i=d"
pknOOS <- read.csv(url1,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
pknOOS$Date <- as.Date(pknOOS$Date)
pknOOS <- pknOOS[, c("Date","Open", "Close")]

url2 <- "http://stooq.com/q/d/l/?s=pgn&i=d"
ltsOOS <- read.csv(url2,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
ltsOOS$Date <- as.Date(ltsOOS$Date)
ltsOOS <- ltsOOS[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20OOS <- read.csv(url1,
                  header = TRUE,
                  sep = ",",
                  dec = ".",
                  stringsAsFactors = F)
wig20OOS$Date <- as.Date(wig20OOS$Date)
wig20OOS <- wig20OOS[, c("Date","Open", "Close")]

pknOOS<-pknOOS[!(pknOOS$Date < "2014-01-01"), ]
pknOOS<-pknOOS[!(pknOOS$Date > "2015-06-30"), ]

ltsOOS<-ltsOOS[!(ltsOOS$Date < "2014-01-01"), ]
ltsOOS<-ltsOOS[!(ltsOOS$Date > "2015-06-30"), ]

wig20OOS<-wig20OOS[!(wig20OOS$Date < "2014-01-01"), ]
wig20OOS<-wig20OOS[!(wig20OOS$Date > "2015-06-30"), ]

pknOOS$Date <- as.Date(pknOOS$Date, format = "%Y-%m-%d")
ltsOOS$Date <- as.Date(ltsOOS$Date, format = "%Y-%m-%d")
wig20OOS$Date <- as.Date(wig20OOS$Date, format = "%Y-%m-%d")

2. Wykres cen zamknięcia:

plot(pknOOS$Date,pknOOS$Close,type="l", ylim=c(4, 20), ylab = "Cena zamknięcia", xlab = "Data", col="lightblue3")
lines(pknOOS$Date, ltsOOS$Close, type="l", col="firebrick3")
title("Ceny zamknięcia Enea i PGNiG (out-of-sample)", cex.main=0.85)

3. Obliczam reszty.

Obliczeń dokonuję na podstawie wektora kointegrującego oszacowanego dla okresu in-sample.

pknOOS$teor=model.coint$coefficients[1]+model.coint$coefficients[2]*ltsOOS$Close
pknOOS$reszty=pknOOS$Close-pknOOS$teor

4. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd. Sprawdzam też, czy dwie kolejne wartości przecięły wartość 0.

pknOOS$czy_sd<-0
pknOOS$czy_zero<-0

for(i in 1:372){
  pknOOS$czy_sd[i] = ifelse(abs(pknOOS$reszty[i])>sd(model.coint$residuals),1,0)
}

for(i in 1:371){
  pknOOS$czy_zero[i+1]=ifelse(pknOOS$reszty[i]*pknOOS$reszty[i+1]>0,0,1)
}

pknOOS$zamk<-0
pknOOS$otw<-0
i=1
k=1

for(k in 1:length(pknOOS$reszty)){
  pknOOS$zamk[k]<-NA
}
for(k in 1:length(pknOOS$reszty)){
  pknOOS$otw[k]<-NA
}

while(i <= 372){
  if(pknOOS$czy_sd[i]==1){
    pknOOS$otw[i]=ifelse(pknOOS$reszty[i]>sd(model.coint$residuals),sd(model.coint$residuals),-sd(model.coint$residuals))
    i<-i+1
    while(i <= 372){
      if(pknOOS$czy_zero[i]==1){
        pknOOS$zamk[i]=0
        break
      }
      i<-i+1
    }
  }
  i<-i+1
}

5. Wykres pozycji dla obu spółek

plot(pknOOS$Date, pknOOS$reszty, type="l", ylab = "Nierównowaga", xlab = "Data")
abline(h=c(sd(model.coint$residuals),-sd(model.coint$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pknOOS$Date, y=pknOOS$otw, col="red", pch=17)
points(x=pknOOS$Date, y=pknOOS$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (out-of-sample)", cex.main=0.85)

6. Obliczam zwroty dla obu spółek.

Następnie zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

i=1
for(i in 2:372){
  pknOOS$zwrot_A[i] = pknOOS$Open[i]/pknOOS$Open[i-1]
}

i=1
for(i in 2:372){
  pknOOS$zwrot_B[i] = ltsOOS$Open[i]/ltsOOS$Open[i-1]
}

pknOOS$otw[pknOOS$otw==sd(model.coint$residuals)]<-1
pknOOS$otw[pknOOS$otw==-sd(model.coint$residuals)]<-1
pknOOS$otw[is.na(pknOOS$otw)]<-0
pknOOS$zamk[pknOOS$zamk==0]<-1
pknOOS$zamk[is.na(pknOOS$zamk)]<-0

7. Tworzę kolumnę pozycja

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1

while(i<=371){
  repeat{
    pknOOS$pozycja[i]=0
    i<-i+1
    if(pknOOS$otw[i-1]!=0||i>371){break}
  }
  while(pknOOS$zamk[i]!=1&&i<=371){
    pknOOS$pozycja[i]=ifelse(pknOOS$reszty[i]>0,1,-1)
    i<-i+1
  }
  pknOOS$pozycja[i]=ifelse(pknOOS$reszty[i-1]>0,1,-1)
  i<-i+1
}

8. Obliczam zwrot z pozycji.

Mówi on, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
pknOOS$zwrot_z_pozycji[1]=0
for(i in 2:372){
  pknOOS$zwrot_z_pozycji[i] = pknOOS$pozycja[i]*(pknOOS$zwrot_B[i]-pknOOS$zwrot_A[i])
}

i=2
for(i in 2:372){
  if(pknOOS$otw[i]==1){
    pknOOS$zwrot_z_pozycji[i+1] = pknOOS$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:372){
  if(pknOOS$zamk[i]==1){
    pknOOS$zwrot_z_pozycji[i] = pknOOS$zwrot_z_pozycji[i]*0.999
  }
}

9. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

pknOOS$portfel[1]=1
i=2
for(i in 2:372){
  pknOOS$portfel[i] = pknOOS$zwrot_z_pozycji[i-1]*pknOOS$portfel[i-1]+pknOOS$portfel[i-1]
}

10. Obliczam dzienne zwroty z indeksu WIG20.

Liczę również wartość portfela dla strategii buy&hold dla WIG20.

i=1
wig20OOS$zwrot<-0
for(i in 2:372){
  wig20OOS$zwrot[i] = wig20OOS$Open[i]/wig20OOS$Open[i-1]
}

wig20OOS$portfel[1]=1
i=2
for(i in 2:372){
  wig20OOS$portfel[i] = wig20OOS$zwrot[i]*wig20OOS$portfel[i-1]
}

11. Wykres equity line.

plot(pknOOS$Date, pknOOS$portfel, type="l", ylim=c(0.9, 1.7), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pknOOS$Date, wig20OOS$portfel, type="l", col="grey")
legend(bty="n", ncol=3,"topleft", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (out-of-sample)", cex.main=.85)

12. Statystyki opisowe dla okresu in-sample

model.coint$zwrot_z_pozycji_na<-model.coint$zwrot_z_pozycji
model.coint$zwrot_z_pozycji_na[model.coint$zwrot_z_pozycji_na==0]<-NA
annual.return<-Return.annualized(na.omit(model.coint$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sd<-sd.annualized(na.omit(model.coint$zwrot_z_pozycji_na), scale = 252)
beta<-lm(model.coint$portfel~wig20$portfel)
drawdown<-maxdrawdown(model.coint$portfel)
wsp.sharpe<-sharpe(model.coint$portfel)
in.sample<-c(annual.return, annual.sd, wsp.sharpe, beta$coefficients[2], drawdown$maxdrawdown)

13. Statystyki opisowe dla okresu out-of-sample.

pknOOS$zwrot_z_pozycji_na<-pknOOS$zwrot_z_pozycji
pknOOS$zwrot_z_pozycji_na[pknOOS$zwrot_z_pozycji_na==0]<-NA
annual.returnOOS<-Return.annualized(na.omit(pknOOS$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sdOOS<-sd.annualized(na.omit(pknOOS$zwrot_z_pozycji_na), scale = 252)
betaOOS<-lm(pknOOS$portfel~wig20OOS$portfel)
drawdownOOS<-maxdrawdown(pknOOS$portfel)
wsp.sharpeOOS<-sharpe(pknOOS$portfel)
wskaźnik<-c("zannualizowana stopa zwrotu", "zannualizowane odchylenie standardowe", "współczynnik Sharpe'a", "współczynnik beta", "maksymalne obsunięcie kapitału")
outof.sample<-c(annual.returnOOS, annual.sdOOS, wsp.sharpeOOS, betaOOS$coefficients[2], drawdownOOS$maxdrawdown)
statystyki<-data.frame(wskaźnik, in.sample, outof.sample)
knitr::kable(statystyki)
wskaźnik in.sample outof.sample
zannualizowana stopa zwrotu -0.0812550 0.3996222
zannualizowane odchylenie standardowe 0.3009941 0.3424174
współczynnik Sharpe’a -0.2630886 0.9613882
współczynnik beta -0.8044628 1.7274482
maksymalne obsunięcie kapitału 0.9014843 0.3561301

Wstępna analiza wykresów equity line pokazuje, że otrzymane wyniki są przeciwne do tych, których można by się spodziewać. W okresie in-sample strategi oparta na kointegracji przynosi straty w porównaniu z okresem out-of-sample.

Wartość zannualizowanej stopy zwrotu w okresie in-sample równa -0,08125499 mówi o tym, że strategia ta przyniosła straty. Takiego wyniku można by się spodziewać dla okresu out-of-sample, dla którego ta wartość wynosi jednak 0.3996222.

Zannualizowane odchylenie standardowe było na podobnym poziomie w obu okresach, w okolicah 30-35%.

Współczynnik Sharpe’a dla okresu in-sample wyniósł -0,2630886. Oznacza to, że zyskowność strategii była mniejsza niż stopa wolna od ryzyka za ten okres. Przeciwne wnioski można wysnuć z wartości współczynnika Sharpe’a dla okresu out-of-sample równego 0.9613882.

Współczynnik beta wyniósł -0.8044628 w okresie in sample, co świadczy o silnej ujemnej korelacji strategii opartej na kointegracji względem buy-and-hold dla indeksu WIG20. Wartośc 1.4876314 dla okresu out-of-sample świadczy z kolei o bardzo silnej, dodatniej korelacji.

Wskaźnika maksymalnego obsunięcia kapitału równy 0.9014843 w okresie in sample, w porównaniu do 0.3561301 świadczy o większym ryzyku związanym ze strategią pair trading w okresie in-sample.

Druga para spółek - PKO BP i ING

Okres in-sample

Wykorzystałem dzienne notowania dla spółek Enea oraz PGNiG. Dane pochodzą z serwisu http://stooq.com.

1. DANE

Importuję dane dotyczące pierwszej spółki oraz dokonuję niezbędnych poprawek: poprawiam format daty oraz usuwam niepotrzebne zmienne:

url1 <- "http://stooq.com/q/d/l/?s=pko&i=d"
pkn2 <- read.csv(url1,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
pkn2$Date <- as.Date(pkn2$Date)
pkn2 <- pkn2[, c("Date","Open", "Close")]

Podobnie postępuję z drugą spółką oraz notowaniami WIG20:

url2 <- "http://stooq.com/q/d/l/?s=ing&i=d"
lts2 <- read.csv(url2,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
lts2$Date <- as.Date(lts2$Date)
lts2 <- lts2[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20 <- read.csv(url1,
                  header = TRUE,
                  sep = ",",
                  dec = ".",
                  stringsAsFactors = F)
wig20$Date <- as.Date(wig20$Date)
wig20 <- wig20[, c("Date","Open", "Close")]

Ograniczam zakres notowań do okresu in-sample oraz przekształcam datę na właściwy format:

pkn2<-pkn2[!(pkn2$Date < "2011-01-01"), ]
pkn2<-pkn2[!(pkn2$Date > "2013-12-31"), ]

lts2<-lts2[!(lts2$Date < "2011-01-01"), ]
lts2<-lts2[!(lts2$Date > "2013-12-31"), ]

wig20<-wig20[!(wig20$Date < "2011-01-01"), ]
wig20<-wig20[!(wig20$Date > "2013-12-31"), ]

pkn2$Date <- as.Date(pkn2$Date, format = "%Y-%m-%d")
lts2$Date <- as.Date(lts2$Date, format = "%Y-%m-%d")
wig20$Date <- as.Date(wig20$Date, format = "%Y-%m-%d")

Wykres cen zamknięcia:

plot(pkn2$Date,pkn2$Close,type="l", ylim=c(20, 115), ylab = "Cena zamknięcia", xlab = "Data", col="lightblue3")
lines(pkn2$Date, lts2$Close, type="l", col="firebrick3")
title("Ceny zamknięcia PKO i ING (in-sample)", cex.main=0.85)

2. Obliczam pierwsze różnice:

lts2$dlts <- diff.xts(lts2$Close)
pkn2$dpkn <- diff.xts(pkn2$Close)

3. Sprawdzam stopień zintegrowania zmiennych

testdf(variable = pkn2$Close,  adf_order = 1)

##   order       adf  p_adf    bgodfrey      p_bg
## 2     0 -2.175589 >10pct 0.006012454 0.9381940
## 3     1 -2.100599 >10pct 0.004085542 0.9490353
testdf(variable = pkn2$dpkn, adf_order = 1)

##   order       adf p_adf    bgodfrey      p_bg
## 2     0 -27.57939 <1pct 0.006860028 0.9339905
## 3     1 -20.50667 <1pct 0.064477147 0.7995547
testdf(variable = lts2$Close,  adf_order = 1)

##   order       adf  p_adf     bgodfrey      p_bg
## 2     0 -1.146521 >10pct 0.4799281158 0.4884549
## 3     1 -1.081743 >10pct 0.0006482768 0.9796870
testdf(variable = lts2$dlts, adf_order = 1)

##   order       adf p_adf    bgodfrey      p_bg
## 2     0 -28.08284 <1pct 0.001181727 0.9725771
## 3     1 -20.90858 <1pct 0.012934160 0.9094531

Statystyka testu Dickeya-Fullera dla spółki PKO BP wynosi -2.175589. Oznacza to brak podstaw do odrzucenia hipotezy zerowej o braku stacjonarności zmiennej losowej. Następnie przeprowadziłem test dla pierwszych róznic cen akcji spółki PKO BP. Wartość statystykli DF wyniosła -27.57939. Oznacza to odrzucenie hipotezy zerowej, że zmienna jest zintegrowana co najmniej w stopniu 2. Wartość statystyki Breuscha-Godfreya równa 0.006860028 oraz p-value wynoszące 0.9339905 świadczą o tym, że realizacja jednokrotnie zróżnicowanej zmiennej jest białym szumem. Świadczy o tym również kształt otrzymanego wykresu.

Wynika z tego, że ceny akcji PKO BP są zintegrowane ~I(1).

Odpowiednie wartości statystyk dla akcji spółki ING wyniosły:

  1. -1.146521
  2. -28.08284
  3. 0.001181727
  4. 0.9725771

Na podstawie powyższych wartości jak również analizy otrzymanego wykresu stwierdzam, że ceny akcji ING są zintegrowane ~I(1)

Obie zmienne są ~I(1) więc sprawdzam, czy są skointegrowane.

4. Szacuję wektor kointegrujący

model.coint2 <- lm(pkn2$Close ~ lts2$Close)
summary(model.coint2)
## 
## Call:
## lm(formula = pkn2$Close ~ lts2$Close)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -4.0456 -1.6682 -0.4145  1.2739  6.3852 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 17.195198   0.617993   27.82   <2e-16 ***
## lts2$Close   0.195358   0.007422   26.32   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.257 on 745 degrees of freedom
## Multiple R-squared:  0.4819, Adjusted R-squared:  0.4812 
## F-statistic: 692.8 on 1 and 745 DF,  p-value: < 2.2e-16

5. Testuję stacjonarność reszt

testdf(variable = model.coint2$residuals, adf_order = 1)

##   order       adf p_adf   bgodfrey       p_bg
## 2     0 -3.267107 <5pct 4.90154762 0.02683264
## 3     1 -2.898313 <5pct 0.04059356 0.84032448

W celu zbadania kointegracji obu zmiennych należy oszacować wektor kointegrujący i przetestować stacjonarność reszt.

Wartość statystyki DF równa -3.267107 oznacza brak podstaw do odrzucenia hipotezy zerowej o stacjonarnosci reszt.

Wartość statystyki BG równa 4.90154762 oraz wartość p-value na poziomie 0.02683264 oznaczają brak podstaw do odrzucenia hipotezy o braku autokorelacji reszt.

Oznacza to, że reszty z oszacowanego modelu są stacjonarne, co świadczy o kointegracji spółek PKO BP i ING.

6. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Skoro ceny obu spółek są skointegrowane, będę mógł zastosować strategię pair-trading. W tym celu muszę wyznaczyć sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd

model.coint2$czy_sd<-0
model.coint2$czy_zero<-0
for(i in 1:747){
  model.coint2$czy_sd[i] = ifelse(abs(model.coint2$residuals[i])>sd(model.coint2$residuals),1,0)
}

Sprawdzam, czy dwie kolejne wartości przecięły wartość 0.

for(i in 1:746){
  model.coint2$czy_zero[i+1]=ifelse(model.coint2$residuals[i]*model.coint2$residuals[i+1]>0,0,1)
}

model.coint2$zamk<-0
model.coint2$otw<-0
i=1
k=1

for(k in 1:length(model.coint2$residuals)){
  model.coint2$zamk[k]<-NA
}
for(k in 1:length(model.coint2$residuals)){
  model.coint2$otw[k]<-NA
}

Uzupełniam kolumny otw oraz zamk sygnałami.

while(i <= 747){
  if(model.coint2$czy_sd[i]==1){
    model.coint2$otw[i]=ifelse(model.coint2$residuals[i]>sd(model.coint2$residuals),sd(model.coint2$residuals),-sd(model.coint2$residuals))
    i<-i+1
    while(i <= 747){
      if(model.coint2$czy_zero[i]==1){
        model.coint2$zamk[i]=0
        break
      }
      i<-i+1
    }
  }
  i<-i+1
}

7. Wykres pozycji dla obu spółek

plot(pkn2$Date, model.coint2$residuals, type="l", ylab = "Nierównowaga", xlab = "Data")
abline(h=c(sd(model.coint2$residuals),-sd(model.coint2$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pkn2$Date, y=model.coint2$otw, col="red", pch=17)
points(x=pkn2$Date, y=model.coint2$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (in-sample)", cex.main=0.85)

Powyższy wykres pokazuje sygnały otwarcia/zamknięcia pozycji w strategii pair-traiding w okresie in-sample. Opierają się one na wykresie wartości nierównowagi krótkookresowej, czyli tak naprawdę na wartościach reszt z oszacowanego modelu. Moment wychylenia się nierównowagi krótkookresowej poza wartość jednego odchylenia standardowego (czerwona linia) stanowi sygnał otwarcia pozycji (czerwony trójkąt). Pozycję zamykam (niebieski kwadrat), kiedy wartość reszty wróci do średniej, czyli wartości 0 (niebieska linia).

8. Obliczam zwroty dla obu spółek.

Zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

i=1
for(i in 2:747){
  model.coint2$zwrot_A[i] = pkn2$Open[i]/pkn2$Open[i-1]
}

i=1
for(i in 2:747){
  model.coint2$zwrot_B[i] = lts2$Open[i]/lts2$Open[i-1]
}

model.coint2$otw[model.coint2$otw==sd(model.coint2$residuals)]<-1
model.coint2$otw[model.coint2$otw==-sd(model.coint2$residuals)]<-1
model.coint2$otw[is.na(model.coint2$otw)]<-0
model.coint2$zamk[model.coint2$zamk==0]<-1
model.coint2$zamk[is.na(model.coint2$zamk)]<-0

9. Tworzę kolumnę pozycja.

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1
while(i<=746){
  repeat{
    model.coint2$pozycja[i]=0
    i<-i+1
    if(model.coint2$otw[i-1]!=0||i>746){break}
  }
  while(model.coint2$zamk[i]!=1&&i<=746){
    model.coint2$pozycja[i]=ifelse(model.coint2$residuals[i]>0,1,-1)
    i<-i+1
  }
  model.coint2$pozycja[i]=ifelse(model.coint2$residuals[i-1]>0,1,-1)
  i<-i+1
}

10. Obliczam zwrot z pozycji.

Kolumna zwrot z pozycji mówi, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
model.coint2$zwrot_z_pozycji[1]=0
for(i in 2:747){
  model.coint2$zwrot_z_pozycji[i] = model.coint2$pozycja[i]*(model.coint2$zwrot_B[i]-model.coint2$zwrot_A[i])
}

i=2
for(i in 2:747){
  if(model.coint2$otw[i]==1){
    model.coint2$zwrot_z_pozycji[i+1] = model.coint2$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:747){
  if(model.coint2$zamk[i]==1){
    model.coint2$zwrot_z_pozycji[i] = model.coint2$zwrot_z_pozycji[i]*0.999
  }
}

11. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

model.coint2$portfel[1]=1
i=2
for(i in 2:747){
  model.coint2$portfel[i] = model.coint2$zwrot_z_pozycji[i-1]*model.coint2$portfel[i-1]+model.coint2$portfel[i-1]
}

12. Obliczam dzienne zwroty z indeksu WIG20

Liczę też zwrot z portfela dla strategii buy&hold dla WIG20.

i=2
wig20$zwrot<-0
for(i in 2:747){
  wig20$zwrot[i] = wig20$Open[i]/wig20$Open[i-1]
}

wig20$portfel[1]=1
i=2
for(i in 2:747){
  wig20$portfel[i] = wig20$zwrot[i]*wig20$portfel[i-1]
}

13. Wykres equity line

plot(pkn2$Date, model.coint2$portfel, type="l", ylim=c(0.7, 1.3), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pkn2$Date, wig20$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "topleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (in-sample)", cex.main=.85)

Praktycznie do połowy roku 2011 widoczna jest wyraźna przewaga zyskowności strategii opartej na kointegracji nad buy&hold na WIG20. Jej rentowność stabilnie utrzymuje się ponad rentownością strategii buy&hold na indeks WIG20

Okres out-of-sample

1. DANE

Importuję dane i dokonuję niezbędnych poprawek na tych samych spółkach w okresie out-of-sample.

url1 <- "http://stooq.com/q/d/l/?s=pko&i=d"
pknOOS2 <- read.csv(url1,
                   header = TRUE,
                   sep = ",",
                   dec = ".",
                   stringsAsFactors = F)
pknOOS2$Date <- as.Date(pknOOS2$Date)
pknOOS2 <- pknOOS2[, c("Date","Open", "Close")]

url2 <- "http://stooq.com/q/d/l/?s=ing&i=d"
ltsOOS2 <- read.csv(url2,
                   header = TRUE,
                   sep = ",",
                   dec = ".",
                   stringsAsFactors = F)
ltsOOS2$Date <- as.Date(ltsOOS2$Date)
ltsOOS2 <- ltsOOS2[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20OOS <- read.csv(url1,
                     header = TRUE,
                     sep = ",",
                     dec = ".",
                     stringsAsFactors = F)
wig20OOS$Date <- as.Date(wig20OOS$Date)
wig20OOS <- wig20OOS[, c("Date","Open", "Close")]

pknOOS2<-pknOOS2[!(pknOOS2$Date < "2014-01-01"), ]
pknOOS2<-pknOOS2[!(pknOOS2$Date > "2015-06-30"), ]

ltsOOS2<-ltsOOS2[!(ltsOOS2$Date < "2014-01-01"), ]
ltsOOS2<-ltsOOS2[!(ltsOOS2$Date > "2015-06-30"), ]

wig20OOS<-wig20OOS[!(wig20OOS$Date < "2014-01-01"), ]
wig20OOS<-wig20OOS[!(wig20OOS$Date > "2015-06-30"), ]

pknOOS2$Date <- as.Date(pknOOS2$Date, format = "%Y-%m-%d")
ltsOOS2$Date <- as.Date(ltsOOS2$Date, format = "%Y-%m-%d")
wig20OOS$Date <- as.Date(wig20OOS$Date, format = "%Y-%m-%d")

2. Wykres cen zamknięcia:

plot(pknOOS2$Date, pknOOS2$Close, type="l", ylab = "Cena zamknięcia", xlab = "Data", ylim=c(30, 160), col="lightblue3")
lines(pknOOS2$Date, ltsOOS2$Close, type="l", col="firebrick3")
title("Ceny zamknięcia PKO i ING (out-of-sample)", cex.main=0.85)

3. Obliczam reszty.

Obliczeń dokonuję na podstawie wektora kointegrującego oszacowanego dla okresu in-sample.

pknOOS2$teor=model.coint2$coefficients[1]+model.coint2$coefficients[2]*ltsOOS2$Close
pknOOS2$reszty=pknOOS2$Close-pknOOS2$teor

4. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd. Sprawdzam też, czy dwie kolejne wartości przecięły wartość 0.

pknOOS2$czy_sd<-0
pknOOS2$czy_zero<-0

for(i in 1:372){
  pknOOS2$czy_sd[i] = ifelse(abs(pknOOS2$reszty[i])>sd(model.coint2$residuals),1,0)
}

for(i in 1:371){
  pknOOS2$czy_zero[i+1]=ifelse(pknOOS2$reszty[i]*pknOOS2$reszty[i+1]>0,0,1)
}

pknOOS2$zamk<-0
pknOOS2$otw<-0
i=1
k=1

for(k in 1:length(pknOOS2$reszty)){
  pknOOS2$zamk[k]<-NA
}
for(k in 1:length(pknOOS2$reszty)){
  pknOOS2$otw[k]<-NA
}

while(i <= 372){
  if(pknOOS2$czy_sd[i]==1){
    pknOOS2$otw[i]=ifelse(pknOOS2$reszty[i]>sd(model.coint2$residuals),sd(model.coint2$residuals),-sd(model.coint2$residuals))
    i<-i+1
    while(i <= 372){
      if(pknOOS2$czy_zero[i]==1){
        pknOOS2$zamk[i]=0
        break
      }
      i<-i+1
    }
  }
  i<-i+1
}

5. Wykres pozycji dla obu spółek

plot(pknOOS2$Date, pknOOS2$reszty, type="l", ylab = "Nierównowaga", xlab = "Data")
abline(h=c(sd(model.coint2$residuals),-sd(model.coint2$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pknOOS2$Date, y=pknOOS2$otw, col="red", pch=17)
points(x=pknOOS2$Date, y=pknOOS2$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (out-of-sample)", cex.main=0.85)

6. Obliczam zwroty dla obu spółek.

Następnie zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

pknOOS2$zwrot_A<-0
pknOOS2$zwrot_B<-0
i=1
for(i in 2:372){
  pknOOS2$zwrot_A[i] = pknOOS2$Open[i]/pknOOS2$Open[i-1]
}

i=1
for(i in 2:372){
  pknOOS2$zwrot_B[i] = ltsOOS2$Open[i]/ltsOOS2$Open[i-1]
}

pknOOS2$otw[pknOOS2$otw==sd(model.coint2$residuals)]<-1
pknOOS2$otw[pknOOS2$otw==-sd(model.coint2$residuals)]<-1
pknOOS2$otw[is.na(pknOOS2$otw)]<-0
pknOOS2$zamk[pknOOS2$zamk==0]<-1
pknOOS2$zamk[is.na(pknOOS2$zamk)]<-0

7. Tworzę kolumnę pozycja

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1

while(i<=371){
  repeat{
    pknOOS2$pozycja[i]=0
    i<-i+1
    if(pknOOS2$otw[i-1]!=0||i>371){break}
  }
  while(pknOOS2$zamk[i]!=1&&i<=371){
    pknOOS2$pozycja[i]=ifelse(pknOOS2$reszty[i]>0,1,-1)
    i<-i+1
  }
  pknOOS2$pozycja[i]=ifelse(pknOOS2$reszty[i-1]>0,1,-1)
  i<-i+1
}

8. Obliczam zwrot z pozycji.

Mówi on, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
pknOOS2$zwrot_z_pozycji[1]=0
for(i in 2:372){
  pknOOS2$zwrot_z_pozycji[i] = pknOOS2$pozycja[i]*(pknOOS2$zwrot_B[i]-pknOOS2$zwrot_A[i])
}

i=2
for(i in 2:372){
  if(pknOOS2$otw[i]==1){
    pknOOS2$zwrot_z_pozycji[i+1] = pknOOS2$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:372){
  if(pknOOS2$zamk[i]==1){
    pknOOS2$zwrot_z_pozycji[i] = pknOOS2$zwrot_z_pozycji[i]*0.999
  }
}

9. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

pknOOS2$portfel[1]=1
i=2
for(i in 2:372){
  pknOOS2$portfel[i] = pknOOS2$zwrot_z_pozycji[i-1]*pknOOS2$portfel[i-1]+pknOOS2$portfel[i-1]
}

10. Obliczam dzienne zwroty z indeksu WIG20.

Liczę również wartość portfela dla strategii buy&hold dla WIG20.

i=1
wig20OOS$zwrot<-0
for(i in 2:372){
  wig20OOS$zwrot[i] = wig20OOS$Open[i]/wig20OOS$Open[i-1]
}

wig20OOS$portfel[1]=1
i=2
for(i in 2:372){
  wig20OOS$portfel[i] = wig20OOS$zwrot[i]*wig20OOS$portfel[i-1]
}

11. Wykres equity line.

plot(pknOOS2$Date, pknOOS2$portfel, type="l", ylim=c(0.75, 1.15), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pknOOS2$Date, wig20OOS$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "bottomleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (out-of-sample)", cex.main=.85)

12. Statystyki opisowe dla okresu in-sample

model.coint2$zwrot_z_pozycji_na<-model.coint2$zwrot_z_pozycji
model.coint2$zwrot_z_pozycji_na[model.coint2$zwrot_z_pozycji_na==0]<-NA
annual.return<-Return.annualized(na.omit(model.coint2$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sd<-sd.annualized(na.omit(model.coint2$zwrot_z_pozycji_na), scale = 252)
beta<-lm(model.coint2$portfel~wig20$portfel)
drawdown<-maxdrawdown(model.coint2$portfel)
wsp.sharpe<-sharpe(model.coint2$portfel)
in.sample<-c(annual.return, annual.sd, wsp.sharpe, beta$coefficients[2], drawdown$maxdrawdown)

13. Statystyki opisowe dla okresu out-of-sample.

pknOOS2$zwrot_z_pozycji_na<-pknOOS2$zwrot_z_pozycji
pknOOS2$zwrot_z_pozycji_na[pknOOS2$zwrot_z_pozycji_na==0]<-NA
annual.returnOOS<-Return.annualized(na.omit(pknOOS2$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sdOOS<-sd.annualized(na.omit(pknOOS2$zwrot_z_pozycji_na), scale = 252)
betaOOS<-lm(pknOOS2$portfel~wig20OOS$portfel)
drawdownOOS<-maxdrawdown(pknOOS2$portfel)
wsp.sharpeOOS<-sharpe(pknOOS2$portfel)
wskaźnik<-c("zannualizowana stopa zwrotu", "zannualizowane odchylenie standardowe", "współczynnik Sharpe'a", "współczynnik beta", "maksymalne obsunięcie kapitału")
outof.sample<-c(annual.returnOOS, annual.sdOOS, wsp.sharpeOOS, betaOOS$coefficients[2], drawdownOOS$maxdrawdown)
statystyki<-data.frame(wskaźnik, in.sample, outof.sample)
knitr::kable(statystyki)
wskaźnik in.sample outof.sample
zannualizowana stopa zwrotu 0.1617381 -0.1358065
zannualizowane odchylenie standardowe 0.3089166 0.3047598
współczynnik Sharpe’a 0.2973430 -0.5332751
współczynnik beta -0.3963540 0.8379209
maksymalne obsunięcie kapitału 0.1717359 0.3911048

Wstępna analiza wykresów equity line pokazuje, że otrzymane wyniki są zgodne z przewidywanymi. W okresie in-sample strategi oparta na kointegracji przynosi zyski w porównaniu z okresem out-of-sample.

Wartość zannualizowanej stopy zwrotu w okresie in-sample równa 0.1595752 mówi o tym, że strategia ta przyniosła zyski. Dla okresu out-of-sample ta wartość wynosi -0.1408709 co oznacza, że przyniosła straty.

Zannualizowane odchylenie standardowe było praktycznie identyczne w obu okresach, w granicach 30%.

Współczynnik Sharpe’a dla okresu in-sample wyniósł 0.2973430. Oznacza to, że zyskowność strategii była większa niż stopa wolna od ryzyka za ten okres. Przeciwne wnioski można wysnuć z wartości współczynnika Sharpe’a dla okresu out-of-sample równego -0.5243374.

Współczynnik beta wyniósł -0.3963540 w okresie in-sample, co świadczy o ujemnej korelacji strategii opartej na kointegracji względem buy-and-hold dla indeksu WIG20. Wartośc 0.8409120 dla okresu out-of-sample świadczy z kolei o silnej, dodatniej korelacji.

Wskaźnika maksymalnego obsunięcia kapitału równy 0.1717359 w okresie in sample, w porównaniu do 0.3911048 świadczy o mniejszym ryzyku związanym ze strategią pair trading w okresie in-sample.

Trzecia para spółek - PGE i Tauron

Okres in-sample

Wykorzystałem dzienne notowania dla spółek PGE oraz Tauron. Dane pochodzą z serwisu http://stooq.com.

1. DANE

Importuję dane dotyczące pierwszej spółki oraz dokonuję niezbędnych poprawek: poprawiam format daty oraz usuwam niepotrzebne zmienne:

url1 <- "http://stooq.com/q/d/l/?s=pge&i=d"
pkn3 <- read.csv(url1,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
pkn3$Date <- as.Date(pkn3$Date)
pkn3 <- pkn3[, c("Date","Open", "Close")]

Podobnie postępuję z drugą spółką oraz notowaniami WIG20:

url2 <- "http://stooq.com/q/d/l/?s=tpe&i=d"
lts3 <- read.csv(url2,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
lts3$Date <- as.Date(lts3$Date)
lts3 <- lts3[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20 <- read.csv(url1,
                  header = TRUE,
                  sep = ",",
                  dec = ".",
                  stringsAsFactors = F)
wig20$Date <- as.Date(wig20$Date)
wig20 <- wig20[, c("Date","Open", "Close")]

Ograniczam zakres notowań do okresu in-sample oraz przekształcam datę na właściwy format:

pkn3<-pkn3[!(pkn3$Date < "2011-01-01"), ]
pkn3<-pkn3[!(pkn3$Date > "2013-12-31"), ]

lts3<-lts3[!(lts3$Date < "2011-01-01"), ]
lts3<-lts3[!(lts3$Date > "2013-12-31"), ]

wig20<-wig20[!(wig20$Date < "2011-01-01"), ]
wig20<-wig20[!(wig20$Date > "2013-12-31"), ]

pkn3$Date <- as.Date(pkn3$Date, format = "%Y-%m-%d")
lts3$Date <- as.Date(lts3$Date, format = "%Y-%m-%d")
wig20$Date <- as.Date(wig20$Date, format = "%Y-%m-%d")

Wykres cen zamknięcia:

plot(pkn3$Date, pkn3$Close, type="l", , ylab = "Cena zamknięcia", xlab = "Data", ylim=c(4,20), col="lightblue3")
lines(pkn3$Date, lts3$Close, col="firebrick3")
title("Ceny zamknięcia PGE i Tauron (in-sample)", cex.main=0.85)

2. Obliczam pierwsze różnice:

lts3$dlts <- diff.xts(lts3$Close)
pkn3$dpkn <- diff.xts(pkn3$Close)

3. Sprawdzam stopień zintegrowania zmiennych

testdf(variable = pkn3$Close,  adf_order = 1)

##   order       adf p_adf     bgodfrey      p_bg
## 2     0 -2.993998 <5pct 0.0460311812 0.8301191
## 3     1 -2.905218 <5pct 0.0003414112 0.9852581
testdf(variable = pkn3$dpkn, adf_order = 1)

##   order       adf p_adf     bgodfrey      p_bg
## 2     0 -27.38483 <1pct 0.0000850319 0.9926426
## 3     1 -20.74569 <1pct 0.0053718818 0.9415729
testdf(variable = lts3$Close,  adf_order = 1)

##   order       adf  p_adf     bgodfrey      p_bg
## 2     0 -2.458015 >10pct 0.5868090970 0.4436557
## 3     1 -2.537448 >10pct 0.0008760996 0.9763869
testdf(variable = lts3$dlts, adf_order = 1)

##   order       adf p_adf     bgodfrey      p_bg
## 2     0 -26.59765 <1pct 0.0007728002 0.9778222
## 3     1 -19.12101 <1pct 0.0009516516 0.9753901

Statystyka testu Dickeya-Fullera dla spółki PGE wynosi -2.993998. Oznacza to brak podstaw do odrzucenia hipotezy zerowej o braku stacjonarności zmiennej losowej. Następnie przeprowadziłem test dla pierwszych róznic cen akcji spółki PGE. Wartość statystykli DF wyniosła -27.38483. Oznacza to odrzucenie hipotezy zerowej, że zmienna jest zintegrowana co najmniej w stopniu 2. Wartość statystyki Breuscha-Godfreya równa 0.0000850319 oraz p-value wynoszące 0.9926426 świadczą o tym, że realizacja jednokrotnie zróżnicowanej zmiennej jest białym szumem. Świadczy o tym również kształt otrzymanego wykresu.

Wynika z tego, że ceny akcji PGE są zintegrowane ~I(1).

Odpowiednie wartości statystyk dla akcji spółki Tauron wyniosły:

  1. -2.458015
  2. -26.59765
  3. 0.0007728002
  4. 0.9778222

Na podstawie powyższych wartości jak również analizy otrzymanego wykresu stwierdzam, że ceny akcji Tauron są zintegrowane ~I(1)

Obie zmienne są ~I(1) więc sprawdzam, czy są skointegrowane.

4. Szacuję wektor kointegrujący

model.coint3 <- lm(pkn3$Close ~ lts3$Close)
summary(model.coint3)
## 
## Call:
## lm(formula = pkn3$Close ~ lts3$Close)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.48528 -0.30203  0.03957  0.45202  1.44961 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  5.92159    0.22569   26.24   <2e-16 ***
## lts3$Close   2.40057    0.05144   46.67   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.6974 on 745 degrees of freedom
## Multiple R-squared:  0.7451, Adjusted R-squared:  0.7448 
## F-statistic:  2178 on 1 and 745 DF,  p-value: < 2.2e-16

5. Testuję stacjonarność reszt

testdf(variable = model.coint3$residuals, adf_order = 1)

##   order       adf p_adf   bgodfrey       p_bg
## 2     0 -4.955945 <1pct 3.10398463 0.07810092
## 3     1 -4.612531 <1pct 0.01447824 0.90422527

W celu zbadania kointegracji obu zmiennych należy oszacować wektor kointegrujący i przetestować stacjonarność reszt.

Wartość statystyki DF równa -4.955945 oznacza brak podstaw do odrzucenia hipotezy zerowej o stacjonarnosci reszt.

Wartość statystyki BG równa 3.10398463 oraz wartość p-value na poziomie 0.07810092 oznaczają brak podstaw do odrzucenia hipotezy o braku autokorelacji reszt.

Oznacza to, że reszty z oszacowanego modelu są stacjonarne, co świadczy o kointegracji spółek PGE i Tauron.

6. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Skoro ceny obu spółek są skointegrowane, będę mógł zastosować strategię pair-trading. W tym celu muszę wyznaczyć sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd

model.coint3$czy_sd<-0
model.coint3$czy_zero<-0
for(i in 1:747){
  model.coint3$czy_sd[i] = ifelse(abs(model.coint3$residuals[i])>sd(model.coint3$residuals),1,0)
}

Sprawdzam, czy dwie kolejne wartości przecięły wartość 0.

for(i in 1:746){
  model.coint3$czy_zero[i+1]=ifelse(model.coint3$residuals[i]*model.coint3$residuals[i+1]>0,0,1)
}

model.coint3$zamk<-0
model.coint3$otw<-0
i=1
k=1

for(k in 1:length(model.coint3$residuals)){
  model.coint3$zamk[k]<-NA
}
for(k in 1:length(model.coint3$residuals)){
  model.coint3$otw[k]<-NA
}

Uzupełniam kolumny otw oraz zamk sygnałami.

while(i <= 747){
  if(model.coint3$czy_sd[i]==1){
    model.coint3$otw[i]=ifelse(model.coint3$residuals[i]>sd(model.coint3$residuals),sd(model.coint3$residuals),-sd(model.coint3$residuals))
    i<-i+1
    while(i <= 747){
      if(model.coint3$czy_zero[i]==1){
        model.coint3$zamk[i]=0
        break
      }
      i<-i+1
    }
  }
  i<-i+1
}

7. Wykres pozycji dla obu spółek

plot(pkn3$Date, model.coint3$residuals, type="l", ylab = "Nierównowaga", xlab = "Data")
abline(h=c(sd(model.coint3$residuals),-sd(model.coint3$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pkn3$Date, y=model.coint3$otw, col="red", pch=17)
points(x=pkn3$Date, y=model.coint3$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (in-sample)", cex.main=0.85)

Powyższy wykres pokazuje sygnały otwarcia/zamknięcia pozycji w strategii pair-traiding w okresie in-sample. Opierają się one na wykresie wartości nierównowagi krótkookresowej, czyli tak naprawdę na wartościach reszt z oszacowanego modelu. Moment wychylenia się nierównowagi krótkookresowej poza wartość jednego odchylenia standardowego (czerwona linia) stanowi sygnał otwarcia pozycji (czerwony trójkąt). Pozycję zamykam (niebieski kwadrat), kiedy wartość reszty wróci do średniej, czyli wartości 0 (niebieska linia).

8. Obliczam zwroty dla obu spółek.

Zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

i=1
for(i in 2:747){
  model.coint3$zwrot_A[i] = pkn3$Open[i]/pkn3$Open[i-1]
}

i=1
for(i in 2:747){
  model.coint3$zwrot_B[i] = lts3$Open[i]/lts3$Open[i-1]
}

model.coint3$otw[model.coint3$otw==sd(model.coint3$residuals)]<-1
model.coint3$otw[model.coint3$otw==-sd(model.coint3$residuals)]<-1
model.coint3$otw[is.na(model.coint3$otw)]<-0
model.coint3$zamk[model.coint3$zamk==0]<-1
model.coint3$zamk[is.na(model.coint3$zamk)]<-0

9. Tworzę kolumnę pozycja.

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1
while(i<=746){
  repeat{
    model.coint3$pozycja[i]=0
    i<-i+1
    if(model.coint3$otw[i-1]!=0||i>746){break}
  }
  while(model.coint3$zamk[i]!=1&&i<=746){
    model.coint3$pozycja[i]=ifelse(model.coint3$residuals[i]>0,1,-1)
    i<-i+1
  }
  model.coint3$pozycja[i]=ifelse(model.coint3$residuals[i-1]>0,1,-1)
  i<-i+1
}

10. Obliczam zwrot z pozycji.

Kolumna zwrot z pozycji mówi, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
model.coint3$zwrot_z_pozycji[1]=0
for(i in 2:747){
  model.coint3$zwrot_z_pozycji[i] = model.coint3$pozycja[i]*(model.coint3$zwrot_B[i]-model.coint3$zwrot_A[i])
}

i=2
for(i in 2:747){
  if(model.coint3$otw[i]==1){
    model.coint3$zwrot_z_pozycji[i+1] = model.coint3$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:747){
  if(model.coint3$zamk[i]==1){
    model.coint3$zwrot_z_pozycji[i] = model.coint3$zwrot_z_pozycji[i]*0.999
  }
}

11. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

model.coint3$portfel[1]=1
i=2
for(i in 2:747){
  model.coint3$portfel[i] = model.coint3$zwrot_z_pozycji[i-1]*model.coint3$portfel[i-1]+model.coint3$portfel[i-1]
}

12. Obliczam dzienne zwroty z indeksu WIG20

Liczę też zwrot z portfela dla strategii buy&hold dla WIG20.

i=2
wig20$zwrot<-0
for(i in 2:747){
  wig20$zwrot[i] = wig20$Open[i]/wig20$Open[i-1]
}

wig20$portfel[1]=1
i=2
for(i in 2:747){
  wig20$portfel[i] = wig20$zwrot[i]*wig20$portfel[i-1]
}

13. Wykres equity line

plot(pkn3$Date, model.coint3$portfel, type="l", ylim=c(0.5, 2), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pkn3$Date, wig20$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "topleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (in-sample)", cex.main=.85)

Praktycznie do połowy roku 2011 widoczna jest wyraźna przewaga zyskowności strategii opartej na kointegracji nad buy&hold na WIG20. Jej rentowność stabilnie rośnie i utrzymuje się ponad rentownością strategii buy&hold na indeks WIG20.

Okres out-of-sample

1. DANE

Importuję dane i dokonuję niezbędnych poprawek na tych samych spółkach w okresie out-of-sample.

url1 <- "http://stooq.com/q/d/l/?s=pge&i=d"
pknOOS3 <- read.csv(url1,
                   header = TRUE,
                   sep = ",",
                   dec = ".",
                   stringsAsFactors = F)
pknOOS3$Date <- as.Date(pknOOS3$Date)
pknOOS3 <- pknOOS3[, c("Date","Open", "Close")]

url2 <- "http://stooq.com/q/d/l/?s=tpe&i=d"
ltsOOS3 <- read.csv(url2,
                   header = TRUE,
                   sep = ",",
                   dec = ".",
                   stringsAsFactors = F)
ltsOOS3$Date <- as.Date(ltsOOS3$Date)
ltsOOS3 <- ltsOOS3[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20OOS <- read.csv(url1,
                     header = TRUE,
                     sep = ",",
                     dec = ".",
                     stringsAsFactors = F)
wig20OOS$Date <- as.Date(wig20OOS$Date)
wig20OOS <- wig20OOS[, c("Date","Open", "Close")]

pknOOS3<-pknOOS3[!(pknOOS3$Date < "2014-01-01"), ]
pknOOS3<-pknOOS3[!(pknOOS3$Date > "2015-06-30"), ]

ltsOOS3<-ltsOOS3[!(ltsOOS3$Date < "2014-01-01"), ]
ltsOOS3<-ltsOOS3[!(ltsOOS3$Date > "2015-06-30"), ]

wig20OOS<-wig20OOS[!(wig20OOS$Date < "2014-01-01"), ]
wig20OOS<-wig20OOS[!(wig20OOS$Date > "2015-06-30"), ]

pknOOS3$Date <- as.Date(pknOOS3$Date, format = "%Y-%m-%d")
ltsOOS3$Date <- as.Date(ltsOOS3$Date, format = "%Y-%m-%d")
wig20OOS$Date <- as.Date(wig20OOS$Date, format = "%Y-%m-%d")

2. Wykres cen zamknięcia:

plot(pknOOS3$Date, pknOOS3$Close, type="l", ylim=c(0,25), ylab = "Cena zamknięcia", xlab = "Data", col="lightblue3")
lines(pknOOS3$Date, ltsOOS3$Close, col="firebrick3")
title("Ceny zamknięcia PGE i Tauron (out-of-sample)", cex.main=0.85)

3. Obliczam reszty.

Obliczeń dokonuję na podstawie wektora kointegrującego oszacowanego dla okresu in-sample.

pknOOS3$teor=model.coint3$coefficients[1]+model.coint3$coefficients[2]*ltsOOS3$Close
pknOOS3$reszty=pknOOS3$Close-pknOOS3$teor

4. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd. Sprawdzam też, czy dwie kolejne wartości przecięły wartość 0.

pknOOS3$czy_sd<-0
pknOOS3$czy_zero<-0

for(i in 1:372){
  pknOOS3$czy_sd[i] = ifelse(abs(pknOOS3$reszty[i])>sd(model.coint3$residuals),1,0)
}

for(i in 1:371){
  pknOOS3$czy_zero[i+1]=ifelse(pknOOS3$reszty[i]*pknOOS3$reszty[i+1]>0,0,1)
}

pknOOS3$zamk<-0
pknOOS3$otw<-0
i=1
k=1

for(k in 1:length(pknOOS3$reszty)){
  pknOOS3$zamk[k]<-NA
}
for(k in 1:length(pknOOS3$reszty)){
  pknOOS3$otw[k]<-NA
}

while(i <= 372){
  if(pknOOS3$czy_sd[i]==1){
    pknOOS3$otw[i]=ifelse(pknOOS3$reszty[i]>sd(model.coint3$residuals),sd(model.coint3$residuals),-sd(model.coint3$residuals))
    i<-i+1
    while(i <= 372){
      if(pknOOS3$czy_zero[i]==1){
        pknOOS3$zamk[i]=0
        break
      }
      i<-i+1
    }
  }
  i<-i+1
}

5. Wykres pozycji dla obu spółek

plot(pknOOS3$Date, pknOOS3$reszty, type="l", ylab = "Nierównowaga", xlab = "Data")
abline(h=c(sd(model.coint3$residuals),-sd(model.coint3$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pknOOS3$Date, y=pknOOS3$otw, col="red", pch=17)
points(x=pknOOS3$Date, y=pknOOS3$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (out-of-sample)", cex.main=0.85)

6. Obliczam zwroty dla obu spółek.

Następnie zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

i=1
for(i in 2:372){
  pknOOS3$zwrot_A[i] = pknOOS3$Open[i]/pknOOS3$Open[i-1]
}

i=1
for(i in 2:372){
  pknOOS3$zwrot_B[i] = ltsOOS3$Open[i]/ltsOOS3$Open[i-1]
}

pknOOS3$otw[pknOOS3$otw==sd(model.coint3$residuals)]<-1
pknOOS3$otw[pknOOS3$otw==-sd(model.coint3$residuals)]<-1
pknOOS3$otw[is.na(pknOOS3$otw)]<-0
pknOOS3$zamk[pknOOS3$zamk==0]<-1
pknOOS3$zamk[is.na(pknOOS3$zamk)]<-0

7. Tworzę kolumnę pozycja

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1

while(i<=371){
  repeat{
    pknOOS3$pozycja[i]=0
    i<-i+1
    if(pknOOS3$otw[i-1]!=0||i>371){break}
  }
  while(pknOOS3$zamk[i]!=1&&i<=371){
    pknOOS3$pozycja[i]=ifelse(pknOOS3$reszty[i]>0,1,-1)
    i<-i+1
  }
  pknOOS3$pozycja[i]=ifelse(pknOOS3$reszty[i-1]>0,1,-1)
  i<-i+1
}

8. Obliczam zwrot z pozycji.

Mówi on, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
pknOOS3$zwrot_z_pozycji[1]=0
for(i in 2:372){
  pknOOS3$zwrot_z_pozycji[i] = pknOOS3$pozycja[i]*(pknOOS3$zwrot_B[i]-pknOOS3$zwrot_A[i])
}

i=2
for(i in 2:372){
  if(pknOOS3$otw[i]==1){
    pknOOS3$zwrot_z_pozycji[i+1] = pknOOS3$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:372){
  if(pknOOS3$zamk[i]==1){
    pknOOS3$zwrot_z_pozycji[i] = pknOOS3$zwrot_z_pozycji[i]*0.999
  }
}

9. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

pknOOS3$portfel[1]=1
i=2
for(i in 2:372){
  pknOOS3$portfel[i] = pknOOS3$zwrot_z_pozycji[i-1]*pknOOS3$portfel[i-1]+pknOOS3$portfel[i-1]
}

10. Obliczam dzienne zwroty z indeksu WIG20.

Liczę również wartość portfela dla strategii buy&hold dla WIG20.

i=1
wig20OOS$zwrot<-0
for(i in 2:372){
  wig20OOS$zwrot[i] = wig20OOS$Open[i]/wig20OOS$Open[i-1]
}

wig20OOS$portfel[1]=1
i=2
for(i in 2:372){
  wig20OOS$portfel[i] = wig20OOS$zwrot[i]*wig20OOS$portfel[i-1]
}

11. Wykres equity line.

plot(pknOOS3$Date, pknOOS3$portfel, type="l", ylim=c(0.8, 1.1), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pknOOS3$Date, wig20OOS$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "bottomleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (out-of-sample)", cex.main=.85)

12. Statystyki opisowe dla okresu in-sample

model.coint3$zwrot_z_pozycji_na<-model.coint3$zwrot_z_pozycji
model.coint3$zwrot_z_pozycji_na[model.coint3$zwrot_z_pozycji_na==0]<-NA
annual.return<-Return.annualized(na.omit(model.coint3$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sd<-sd.annualized(na.omit(model.coint3$zwrot_z_pozycji_na), scale = 252)
beta<-lm(model.coint3$portfel~wig20$portfel)
drawdown<-maxdrawdown(model.coint3$portfel)
wsp.sharpe<-sharpe(model.coint3$portfel)
in.sample<-c(annual.return, annual.sd, wsp.sharpe, beta$coefficients[2], drawdown$maxdrawdown)

13. Statystyki opisowe dla okresu out-of-sample.

pknOOS3$zwrot_z_pozycji_na<-pknOOS3$zwrot_z_pozycji
pknOOS3$zwrot_z_pozycji_na[pknOOS3$zwrot_z_pozycji_na==0]<-NA
annual.returnOOS<-Return.annualized(na.omit(pknOOS3$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sdOOS<-sd.annualized(na.omit(pknOOS3$zwrot_z_pozycji_na), scale = 252)
betaOOS<-lm(pknOOS3$portfel~wig20OOS$portfel)
drawdownOOS<-maxdrawdown(pknOOS3$portfel)
wsp.sharpeOOS<-sharpe(pknOOS3$portfel)
wskaźnik<-c("zannualizowana stopa zwrotu", "zannualizowane odchylenie standardowe", "współczynnik Sharpe'a", "współczynnik beta", "maksymalne obsunięcie kapitału")
outof.sample<-c(annual.returnOOS, annual.sdOOS, wsp.sharpeOOS, betaOOS$coefficients[2], drawdownOOS$maxdrawdown)
statystyki<-data.frame(wskaźnik, in.sample, outof.sample)
knitr::kable(statystyki)
wskaźnik in.sample outof.sample
zannualizowana stopa zwrotu 0.4453060 -0.0512173
zannualizowane odchylenie standardowe 0.2854798 0.2443171
współczynnik Sharpe’a 0.9874263 -0.3061049
współczynnik beta -1.6496693 -0.1901901
maksymalne obsunięcie kapitału 0.3150918 0.2795250

Wstępna analiza wykresów equity line pokazuje, że otrzymane wyniki są zgodne z przewidywanymi. W okresie in-sample strategi oparta na kointegracji przynosi zyski w porównaniu z okresem out-of-sample.

Wartość zannualizowanej stopy zwrotu w okresie in-sample równa 0.4453060 mówi o tym, że strategia ta przyniosła zyski. Dla okresu out-of-sample ta wartość wynosi -0.0512173 co oznacza, że przyniosła straty.

Zannualizowane odchylenie standardowe dla okresu in-sample wyniosło 28.55%, dla out-of-sample 24.43%

Współczynnik Sharpe’a dla okresu in-sample wyniósł 0.9874263. Oznacza to, że zyskowność strategii była większa niż stopa wolna od ryzyka za ten okres. Przeciwne wnioski można wysnuć z wartości współczynnika Sharpe’a dla okresu out-of-sample równego -0.3061049.

Współczynnik beta wyniósł -1.6496693 w okresie in-sample, co świadczy o silnej ujemnej korelacji strategii opartej na kointegracji względem buy-and-hold dla indeksu WIG20. Wartośc -0.1901901 dla okresu out-of-sample świadczy o słabej ujemnej korelacji względem indeksu WIG20.

Słabe wyniki strategii w okresie out-of-sample mogą być spowodowane bardzo małą liczbą wykonanych transakcji, co często ogranicza zyskowność tej strategii.

Czwarta para spółek - PKN Orlen i Lotos

Okres in-sample

Wykorzystałem dzienne notowania dla spółek PKN Orlen oraz Lotos. Dane pochodzą z serwisu http://stooq.com.

1. DANE

Importuję dane dotyczące pierwszej spółki oraz dokonuję niezbędnych poprawek: poprawiam format daty oraz usuwam niepotrzebne zmienne:

url1 <- "http://stooq.com/q/d/l/?s=pkn&i=d"
pkn4 <- read.csv(url1,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
pkn4$Date <- as.Date(pkn4$Date)
pkn4 <- pkn4[, c("Date","Open", "Close")]

Podobnie postępuję z drugą spółką oraz notowaniami WIG20:

url2 <- "http://stooq.com/q/d/l/?s=lts&i=d"
lts4 <- read.csv(url2,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
lts4$Date <- as.Date(lts4$Date)
lts4 <- lts4[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20 <- read.csv(url1,
                  header = TRUE,
                  sep = ",",
                  dec = ".",
                  stringsAsFactors = F)
wig20$Date <- as.Date(wig20$Date)
wig20 <- wig20[, c("Date","Open", "Close")]

Ograniczam zakres notowań do okresu in-sample oraz przekształcam datę na właściwy format:

pkn4<-pkn4[!(pkn4$Date < "2011-01-01"), ]
pkn4<-pkn4[!(pkn4$Date > "2013-12-31"), ]

lts4<-lts4[!(lts4$Date < "2011-01-01"), ]
lts4<-lts4[!(lts4$Date > "2013-12-31"), ]

wig20<-wig20[!(wig20$Date < "2011-01-01"), ]
wig20<-wig20[!(wig20$Date > "2013-12-31"), ]

pkn4$Date <- as.Date(pkn4$Date, format = "%Y-%m-%d")
lts4$Date <- as.Date(lts4$Date, format = "%Y-%m-%d")
wig20$Date <- as.Date(wig20$Date, format = "%Y-%m-%d")

Wykres cen zamknięcia:

plot(pkn4$Date, pkn4$Close, type="l", ylim=c(15,65), ylab = "Cena zamknięcia", xlab = "Data", col="lightblue3")
lines(pkn4$Date, lts4$Close, col="firebrick3")
title("Ceny zamknięcia PKN Orlen i Lotos (in-sample)", cex.main=0.85)

2. Obliczam pierwsze różnice:

lts4$dlts <- diff.xts(lts4$Close)
pkn4$dpkn <- diff.xts(pkn4$Close)

3. Sprawdzam stopień zintegrowania zmiennych

testdf(variable = pkn4$Close,  adf_order = 1)

##   order       adf  p_adf     bgodfrey      p_bg
## 2     0 -1.909782 >10pct 0.5323054979 0.4656392
## 3     1 -1.929410 >10pct 0.0003448212 0.9851846
testdf(variable = pkn4$dpkn, adf_order = 1)

##   order       adf p_adf     bgodfrey      p_bg
## 2     0 -26.71328 <1pct 0.0004424368 0.9832184
## 3     1 -19.66189 <1pct 0.0002246286 0.9880421
testdf(variable = lts4$Close,  adf_order = 1)

##   order       adf  p_adf    bgodfrey      p_bg
## 2     0 -1.421652 >10pct 2.683260172 0.1014079
## 3     1 -1.500401 >10pct 0.001989564 0.9644225
testdf(variable = lts4$dlts, adf_order = 1)

##   order       adf p_adf    bgodfrey      p_bg
## 2     0 -25.73044 <1pct 0.001520219 0.9688984
## 3     1 -18.34315 <1pct 0.001597862 0.9681144

Statystyka testu Dickeya-Fullera dla spółki PKN Orlen wynosi -1.909782. Oznacza to brak podstaw do odrzucenia hipotezy zerowej o braku stacjonarności zmiennej losowej. Następnie przeprowadziłem test dla pierwszych róznic cen akcji spółki PKN Orlen. Wartość statystykli DF wyniosła -26.71328. Oznacza to odrzucenie hipotezy zerowej, że zmienna jest zintegrowana co najmniej w stopniu 2. Wartość statystyki Breuscha-Godfreya równa 0.0004424368 oraz p-value wynoszące 0.9832184 świadczą o tym, że realizacja jednokrotnie zróżnicowanej zmiennej jest białym szumem. Świadczy o tym również kształt otrzymanego wykresu.

Wynika z tego, że ceny akcji PKN Orlen są zintegrowane ~I(1).

Odpowiednie wartości statystyk dla akcji spółki Lotos wyniosły:

  1. -1.421652
  2. -25.73044
  3. 0.001520219
  4. 0.9688984

Na podstawie powyższych wartości jak również analizy otrzymanego wykresu stwierdzam, że ceny akcji Lotos są zintegrowane ~I(1)

Obie zmienne są ~I(1) więc sprawdzam, czy są skointegrowane.

4. Szacuję wektor kointegrujący

model.coint4 <- lm(pkn4$Close ~ lts4$Close)
summary(model.coint4)
## 
## Call:
## lm(formula = pkn4$Close ~ lts4$Close)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -6.8060 -1.2375  0.0975  1.2459  5.4052 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 13.30401    0.34534   38.52   <2e-16 ***
## lts4$Close   0.87013    0.01085   80.18   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 1.917 on 745 degrees of freedom
## Multiple R-squared:  0.8962, Adjusted R-squared:  0.896 
## F-statistic:  6429 on 1 and 745 DF,  p-value: < 2.2e-16

5. Testuję stacjonarność reszt

testdf(variable = model.coint4$residuals, adf_order = 1)

##   order       adf p_adf     bgodfrey      p_bg
## 2     0 -4.576841 <1pct 0.0025193495 0.9599685
## 3     1 -4.477310 <1pct 0.0000684903 0.9933969

W celu zbadania kointegracji obu zmiennych należy oszacować wektor kointegrujący i przetestować stacjonarność reszt.

Wartość statystyki DF równa -4.576841 oznacza brak podstaw do odrzucenia hipotezy zerowej o stacjonarnosci reszt.

Wartość statystyki BG równa 0.0025193495 oraz wartość p-value na poziomie 0.9599685 oznaczają brak podstaw do odrzucenia hipotezy o braku autokorelacji reszt.

Oznacza to, że reszty z oszacowanego modelu są stacjonarne, co świadczy o kointegracji spółek PKN Orlen i Lotos.

6. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Skoro ceny obu spółek są skointegrowane, będę mógł zastosować strategię pair-trading. W tym celu muszę wyznaczyć sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd

model.coint4$czy_sd<-0
model.coint4$czy_zero<-0
for(i in 1:747){
  model.coint4$czy_sd[i] = ifelse(abs(model.coint4$residuals[i])>sd(model.coint4$residuals),1,0)
}

Sprawdzam, czy dwie kolejne wartości przecięły wartość 0.

for(i in 1:746){
  model.coint4$czy_zero[i+1]=ifelse(model.coint4$residuals[i]*model.coint4$residuals[i+1]>0,0,1)
}

model.coint4$zamk<-0
model.coint4$otw<-0
i=1
k=1

for(k in 1:length(model.coint4$residuals)){
  model.coint4$zamk[k]<-NA
}
for(k in 1:length(model.coint4$residuals)){
  model.coint4$otw[k]<-NA
}

Uzupełniam kolumny otw oraz zamk sygnałami.

while(i <= 747){
  if(model.coint4$czy_sd[i]==1){
    model.coint4$otw[i]=ifelse(model.coint4$residuals[i]>sd(model.coint4$residuals),sd(model.coint4$residuals),-sd(model.coint4$residuals))
    i<-i+1
    while(i <= 747){
      if(model.coint4$czy_zero[i]==1){
        model.coint4$zamk[i]=0
        break
      }
      i<-i+1
    }
  }
  i<-i+1
}

7. Wykres pozycji dla obu spółek

plot(pkn4$Date, model.coint4$residuals, type="l", ylab = "Nierównowaga", xlab = "Data")
abline(h=c(sd(model.coint4$residuals),-sd(model.coint4$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pkn4$Date, y=model.coint4$otw, col="red", pch=17)
points(x=pkn4$Date, y=model.coint4$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (in-sample)", cex.main=0.85)

Powyższy wykres pokazuje sygnały otwarcia/zamknięcia pozycji w strategii pair-traiding w okresie in-sample. Opierają się one na wykresie wartości nierównowagi krótkookresowej, czyli tak naprawdę na wartościach reszt z oszacowanego modelu. Moment wychylenia się nierównowagi krótkookresowej poza wartość jednego odchylenia standardowego (czerwona linia) stanowi sygnał otwarcia pozycji (czerwony trójkąt). Pozycję zamykam (niebieski kwadrat), kiedy wartość reszty wróci do średniej, czyli wartości 0 (niebieska linia).

8. Obliczam zwroty dla obu spółek.

Zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

i=1
for(i in 2:747){
  model.coint4$zwrot_A[i] = pkn4$Open[i]/pkn4$Open[i-1]
}

i=1
for(i in 2:747){
  model.coint4$zwrot_B[i] = lts4$Open[i]/lts4$Open[i-1]
}

model.coint4$otw[model.coint4$otw==sd(model.coint4$residuals)]<-1
model.coint4$otw[model.coint4$otw==-sd(model.coint4$residuals)]<-1
model.coint4$otw[is.na(model.coint4$otw)]<-0
model.coint4$zamk[model.coint4$zamk==0]<-1
model.coint4$zamk[is.na(model.coint4$zamk)]<-0

9. Tworzę kolumnę pozycja.

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1
while(i<=746){
  repeat{
    model.coint4$pozycja[i]=0
    i<-i+1
    if(model.coint4$otw[i-1]!=0||i>746){break}
  }
  while(model.coint4$zamk[i]!=1&&i<=746){
    model.coint4$pozycja[i]=ifelse(model.coint4$residuals[i]>0,1,-1)
    i<-i+1
  }
  model.coint4$pozycja[i]=ifelse(model.coint4$residuals[i-1]>0,1,-1)
  i<-i+1
}

10. Obliczam zwrot z pozycji.

Kolumna zwrot z pozycji mówi, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
model.coint4$zwrot_z_pozycji[1]=0
for(i in 2:747){
  model.coint4$zwrot_z_pozycji[i] = model.coint4$pozycja[i]*(model.coint4$zwrot_B[i]-model.coint4$zwrot_A[i])
}

i=2
for(i in 2:747){
  if(model.coint4$otw[i]==1){
    model.coint4$zwrot_z_pozycji[i+1] = model.coint4$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:747){
  if(model.coint4$zamk[i]==1){
    model.coint4$zwrot_z_pozycji[i] = model.coint4$zwrot_z_pozycji[i]*0.999
  }
}

11. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

model.coint4$portfel[1]=1
i=2
for(i in 2:747){
  model.coint4$portfel[i] = model.coint4$zwrot_z_pozycji[i-1]*model.coint4$portfel[i-1]+model.coint4$portfel[i-1]
}

12. Obliczam dzienne zwroty z indeksu WIG20

Liczę też zwrot z portfela dla strategii buy&hold dla WIG20.

i=2
wig20$zwrot<-0
for(i in 2:747){
  wig20$zwrot[i] = wig20$Open[i]/wig20$Open[i-1]
}

wig20$portfel[1]=1
i=2
for(i in 2:747){
  wig20$portfel[i] = wig20$zwrot[i]*wig20$portfel[i-1]
}

13. Wykres equity line

plot(pkn4$Date, model.coint4$portfel, type="l", ylim=c(0.5, 2), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pkn4$Date, wig20$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "topleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (in-sample)", cex.main=.85)

Praktycznie do połowy roku 2011 widoczna jest wyraźna przewaga zyskowności strategii opartej na kointegracji nad buy&hold na WIG20. Jej rentowność stabilnie rośnie i utrzymuje się ponad rentownością strategii buy&hold na indeks WIG20.

Okres out-of-sample

1. DANE

Importuję dane i dokonuję niezbędnych poprawek na tych samych spółkach w okresie out-of-sample.

url1 <- "http://stooq.com/q/d/l/?s=pkn&i=d"
pknOOS4 <- read.csv(url1,
                   header = TRUE,
                   sep = ",",
                   dec = ".",
                   stringsAsFactors = F)
pknOOS4$Date <- as.Date(pknOOS4$Date)
pknOOS4 <- pknOOS4[, c("Date","Open", "Close")]

url2 <- "http://stooq.com/q/d/l/?s=lts&i=d"
ltsOOS4 <- read.csv(url2,
                   header = TRUE,
                   sep = ",",
                   dec = ".",
                   stringsAsFactors = F)
ltsOOS4$Date <- as.Date(ltsOOS4$Date)
ltsOOS4 <- ltsOOS4[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20OOS <- read.csv(url1,
                     header = TRUE,
                     sep = ",",
                     dec = ".",
                     stringsAsFactors = F)
wig20OOS$Date <- as.Date(wig20OOS$Date)
wig20OOS <- wig20OOS[, c("Date","Open", "Close")]

pknOOS4<-pknOOS4[!(pknOOS4$Date < "2014-01-01"), ]
pknOOS4<-pknOOS4[!(pknOOS4$Date > "2015-06-30"), ]

ltsOOS4<-ltsOOS4[!(ltsOOS4$Date < "2014-01-01"), ]
ltsOOS4<-ltsOOS4[!(ltsOOS4$Date > "2015-06-30"), ]

wig20OOS<-wig20OOS[!(wig20OOS$Date < "2014-01-01"), ]
wig20OOS<-wig20OOS[!(wig20OOS$Date > "2015-06-30"), ]

pknOOS4$Date <- as.Date(pknOOS4$Date, format = "%Y-%m-%d")
ltsOOS4$Date <- as.Date(ltsOOS4$Date, format = "%Y-%m-%d")
wig20OOS$Date <- as.Date(wig20OOS$Date, format = "%Y-%m-%d")

2. Wykres cen zamknięcia:

plot(pknOOS4$Date, pknOOS4$Close, type="l", ylim=c(20,80), ylab = "Cena zamknięcia", xlab = "Data", col="lightblue3")
lines(pknOOS4$Date, ltsOOS4$Close, col="firebrick3")
title("Ceny zamknięcia PKN Orlen i LOTOS (out-of-sample)", cex.main=0.85)

3. Obliczam reszty.

Obliczeń dokonuję na podstawie wektora kointegrującego oszacowanego dla okresu in-sample.

pknOOS4$teor=model.coint4$coefficients[1]+model.coint4$coefficients[2]*ltsOOS4$Close
pknOOS4$reszty=pknOOS4$Close-pknOOS4$teor

4. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd. Sprawdzam też, czy dwie kolejne wartości przecięły wartość 0.

pknOOS4$czy_sd<-0
pknOOS4$czy_zero<-0

for(i in 1:372){
  pknOOS4$czy_sd[i] = ifelse(abs(pknOOS4$reszty[i])>sd(model.coint4$residuals),1,0)
}

for(i in 1:371){
  pknOOS4$czy_zero[i+1]=ifelse(pknOOS4$reszty[i]*pknOOS4$reszty[i+1]>0,0,1)
}

pknOOS4$zamk<-0
pknOOS4$otw<-0
i=1
k=1

for(k in 1:length(pknOOS4$reszty)){
  pknOOS4$zamk[k]<-NA
}
for(k in 1:length(pknOOS4$reszty)){
  pknOOS4$otw[k]<-NA
}

while(i <= 372){
  if(pknOOS4$czy_sd[i]==1){
    pknOOS4$otw[i]=ifelse(pknOOS4$reszty[i]>sd(model.coint4$residuals),sd(model.coint4$residuals),-sd(model.coint4$residuals))
    i<-i+1
    while(i <= 372){
      if(pknOOS4$czy_zero[i]==1){
        pknOOS4$zamk[i]=0
        break
      }
      i<-i+1
    }
  }
  i<-i+1
}

5. Wykres pozycji dla obu spółek

plot(pknOOS4$Date, pknOOS4$reszty, type="l", ylab = "Nierównowaga", xlab = "Data")
abline(h=c(sd(model.coint4$residuals),-sd(model.coint4$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pknOOS4$Date, y=pknOOS4$otw, col="red", pch=17)
points(x=pknOOS4$Date, y=pknOOS4$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (out-of-sample)", cex.main=0.85)

6. Obliczam zwroty dla obu spółek.

Następnie zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

i=1
for(i in 2:372){
  pknOOS4$zwrot_A[i] = pknOOS4$Open[i]/pknOOS4$Open[i-1]
}

i=1
for(i in 2:372){
  pknOOS4$zwrot_B[i] = ltsOOS4$Open[i]/ltsOOS4$Open[i-1]
}

pknOOS4$otw[pknOOS4$otw==sd(model.coint4$residuals)]<-1
pknOOS4$otw[pknOOS4$otw==-sd(model.coint4$residuals)]<-1
pknOOS4$otw[is.na(pknOOS4$otw)]<-0
pknOOS4$zamk[pknOOS4$zamk==0]<-1
pknOOS4$zamk[is.na(pknOOS4$zamk)]<-0

7. Tworzę kolumnę pozycja

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1
while(i<=371){
  repeat{
    pknOOS4$pozycja[i]=0
    i<-i+1
    if(pknOOS4$otw[i-1]!=0||i>371){break}
  }
  while(pknOOS4$zamk[i]!=1&&i<=371){
    pknOOS4$pozycja[i]=ifelse(pknOOS4$reszty[i]>0,1,-1)
    i<-i+1
  }
  pknOOS4$pozycja[i]=ifelse(pknOOS4$reszty[i-1]>0,1,-1)
  i<-i+1
}

8. Obliczam zwrot z pozycji.

Mówi on, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
pknOOS4$zwrot_z_pozycji[1]=0
for(i in 2:372){
  pknOOS4$zwrot_z_pozycji[i] = pknOOS4$pozycja[i]*(pknOOS4$zwrot_B[i]-pknOOS4$zwrot_A[i])
}

i=2
for(i in 2:372){
  if(pknOOS4$otw[i]==1){
    pknOOS4$zwrot_z_pozycji[i+1] = pknOOS4$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:372){
  if(pknOOS4$zamk[i]==1){
    pknOOS4$zwrot_z_pozycji[i] = pknOOS4$zwrot_z_pozycji[i]*0.999
  }
}

9. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

pknOOS4$portfel[1]=1
i=2
for(i in 2:372){
  pknOOS4$portfel[i] = pknOOS4$zwrot_z_pozycji[i-1]*pknOOS4$portfel[i-1]+pknOOS4$portfel[i-1]
}

10. Obliczam dzienne zwroty z indeksu WIG20.

Liczę również wartość portfela dla strategii buy&hold dla WIG20.

i=1
wig20OOS$zwrot<-0
for(i in 2:372){
  wig20OOS$zwrot[i] = wig20OOS$Open[i]/wig20OOS$Open[i-1]
}

wig20OOS$portfel[1]=1
i=2
for(i in 2:372){
  wig20OOS$portfel[i] = wig20OOS$zwrot[i]*wig20OOS$portfel[i-1]
}

11. Wykres equity line.

plot(pknOOS4$Date, pknOOS4$portfel, type="l", ylim=c(0.6, 1.2), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pknOOS4$Date, wig20OOS$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "topleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (out-of-sample)", cex.main=.85)

12. Statystyki opisowe dla okresu in-sample

model.coint4$zwrot_z_pozycji_na<-model.coint4$zwrot_z_pozycji
model.coint4$zwrot_z_pozycji_na[model.coint4$zwrot_z_pozycji_na==0]<-NA
annual.return<-Return.annualized(na.omit(model.coint4$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sd<-sd.annualized(na.omit(model.coint4$zwrot_z_pozycji_na), scale = 252)
beta<-lm(model.coint4$portfel~wig20$portfel)
drawdown<-maxdrawdown(model.coint4$portfel)
wsp.sharpe<-sharpe(model.coint4$portfel)
in.sample<-c(annual.return, annual.sd, wsp.sharpe, beta$coefficients[2], drawdown$maxdrawdown)

13. Statystyki opisowe dla okresu out-of-sample.

pknOOS4$zwrot_z_pozycji_na<-pknOOS4$zwrot_z_pozycji
pknOOS4$zwrot_z_pozycji_na[pknOOS4$zwrot_z_pozycji_na==0]<-NA
annual.returnOOS<-Return.annualized(na.omit(pknOOS4$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sdOOS<-sd.annualized(na.omit(pknOOS4$zwrot_z_pozycji_na), scale = 252)
betaOOS<-lm(pknOOS4$portfel~wig20OOS$portfel)
drawdownOOS<-maxdrawdown(pknOOS4$portfel)
wsp.sharpeOOS<-sharpe(pknOOS4$portfel)
wskaźnik<-c("zannualizowana stopa zwrotu", "zannualizowane odchylenie standardowe", "współczynnik Sharpe'a", "współczynnik beta", "maksymalne obsunięcie kapitału")
outof.sample<-c(annual.returnOOS, annual.sdOOS, wsp.sharpeOOS, betaOOS$coefficients[2], drawdownOOS$maxdrawdown)
statystyki<-data.frame(wskaźnik, in.sample, outof.sample)
knitr::kable(statystyki)
wskaźnik in.sample outof.sample
zannualizowana stopa zwrotu 0.3776479 -0.2910137
zannualizowane odchylenie standardowe 0.3088333 0.3241910
współczynnik Sharpe’a 0.7353650 -0.9275611
współczynnik beta -1.8551552 1.2184185
maksymalne obsunięcie kapitału 0.2020833 0.5020408

Wstępna analiza wykresów equity line pokazuje, że otrzymane wyniki są zgodne z przewidywanymi. W okresie in-sample strategi oparta na kointegracji przynosi zyski w porównaniu z okresem out-of-sample.

Wartość zannualizowanej stopy zwrotu w okresie in-sample równa 0.3776479 mówi o tym, że strategia ta przyniosła zyski. Dla okresu out-of-sample ta wartość wynosi -0.2910137 co oznacza, że przyniosła straty.

Zannualizowane odchylenie standardowe dla okresu in-sample wyniosło 30.89%, dla out-of-sample 32.42%

Współczynnik Sharpe’a dla okresu in-sample wyniósł 0.7353650. Oznacza to, że zyskowność strategii była większa niż stopa wolna od ryzyka za ten okres. Przeciwne wnioski można wysnuć z wartości współczynnika Sharpe’a dla okresu out-of-sample równego -0.9275611.

Współczynnik beta wyniósł -1.8551552 w okresie in-sample, co świadczy o silnej ujemnej korelacji strategii opartej na kointegracji względem buy-and-hold dla indeksu WIG20. Wartośc 1.2184185 dla okresu out-of-sample świadczy o silnej dodatniej korelacji względem indeksu WIG20.

Wskaźnika maksymalnego obsunięcia kapitału równy 0.2020833 w okresie in sample, w porównaniu do 0.5020408 świadczy o mniejszym ryzyku związanym ze strategią pair trading w okresie in-sample.

Piąta para spółek - ING i Pekao SA

Okres in-sample

Wykorzystałem dzienne notowania dla spółek ING i Pekao SA. Dane pochodzą z serwisu http://stooq.com.

1. DANE

Importuję dane dotyczące pierwszej spółki oraz dokonuję niezbędnych poprawek: poprawiam format daty oraz usuwam niepotrzebne zmienne:

url1 <- "http://stooq.com/q/d/l/?s=ing&i=d"
pkn5 <- read.csv(url1,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
pkn5$Date <- as.Date(pkn5$Date)
pkn5 <- pkn5[, c("Date","Open", "Close")]

Podobnie postępuję z drugą spółką oraz notowaniami WIG20:

url2 <- "http://stooq.com/q/d/l/?s=peo&i=d"
lts5 <- read.csv(url2,
                header = TRUE,
                sep = ",",
                dec = ".",
                stringsAsFactors = F)
lts5$Date <- as.Date(lts5$Date)
lts5 <- lts5[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20 <- read.csv(url1,
                  header = TRUE,
                  sep = ",",
                  dec = ".",
                  stringsAsFactors = F)
wig20$Date <- as.Date(wig20$Date)
wig20 <- wig20[, c("Date","Open", "Close")]

Ograniczam zakres notowań do okresu in-sample oraz przekształcam datę na właściwy format:

pkn5<-pkn5[!(pkn5$Date < "2011-01-01"), ]
pkn5<-pkn5[!(pkn5$Date > "2013-12-31"), ]

lts5<-lts5[!(lts5$Date < "2011-01-01"), ]
lts5<-lts5[!(lts5$Date > "2013-12-31"), ]

wig20<-wig20[!(wig20$Date < "2011-01-01"), ]
wig20<-wig20[!(wig20$Date > "2013-12-31"), ]

pkn5$Date <- as.Date(pkn5$Date, format = "%Y-%m-%d")
lts5$Date <- as.Date(lts5$Date, format = "%Y-%m-%d")
wig20$Date <- as.Date(wig20$Date, format = "%Y-%m-%d")

Wykres cen zamknięcia:

plot(pkn5$Date, pkn5$Close, type="l", ylim=c(60,180), ylab = "Cena zamknięcia", xlab = "Data", col="lightblue3")
lines(pkn5$Date, lts5$Close, col="firebrick3")
title("Ceny zamknięcia ING i Pekao (in-sample)", cex.main=0.85)

2. Obliczam pierwsze różnice:

lts5$dlts <- diff.xts(lts5$Close)
pkn5$dpkn <- diff.xts(pkn5$Close)

3. Sprawdzam stopień zintegrowania zmiennych

testdf(variable = pkn5$Close,  adf_order = 1)

##   order       adf  p_adf     bgodfrey      p_bg
## 2     0 -1.146521 >10pct 0.4799281158 0.4884549
## 3     1 -1.081743 >10pct 0.0006482768 0.9796870
testdf(variable = pkn5$dpkn, adf_order = 1)

##   order       adf p_adf    bgodfrey      p_bg
## 2     0 -28.08284 <1pct 0.001181727 0.9725771
## 3     1 -20.90858 <1pct 0.012934160 0.9094531
testdf(variable = lts5$Close,  adf_order = 1)

##   order      adf  p_adf    bgodfrey      p_bg
## 2     0 -1.86097 >10pct 1.587038290 0.2077498
## 3     1 -1.73939 >10pct 0.006782032 0.9343659
testdf(variable = lts5$dlts, adf_order = 1)

##   order       adf p_adf   bgodfrey     p_bg
## 2     0 -28.73485 <1pct 0.01091278 0.916801
## 3     1 -21.43432 <1pct 0.13134364 0.717043

Statystyka testu Dickeya-Fullera dla spółki ING wynosi -1.146521. Oznacza to brak podstaw do odrzucenia hipotezy zerowej o braku stacjonarności zmiennej losowej. Następnie przeprowadziłem test dla pierwszych róznic cen akcji spółki ING. Wartość statystykli DF wyniosła -28.08284. Oznacza to odrzucenie hipotezy zerowej, że zmienna jest zintegrowana co najmniej w stopniu 2. Wartość statystyki Breuscha-Godfreya równa 0.001181727 oraz p-value wynoszące 0.9725771 świadczą o tym, że realizacja jednokrotnie zróżnicowanej zmiennej jest białym szumem. Świadczy o tym również kształt otrzymanego wykresu.

Wynika z tego, że ceny akcji ING są zintegrowane ~I(1).

Odpowiednie wartości statystyk dla akcji spółki Pekao SA wyniosły:

  1. -1.86097
  2. -28.73485
  3. 0.01091278
  4. 0.916801

Na podstawie powyższych wartości jak również analizy otrzymanego wykresu stwierdzam, że ceny akcji Pekao SA są zintegrowane ~I(1)

Obie zmienne są ~I(1) więc sprawdzam, czy są skointegrowane.

4. Szacuję wektor kointegrujący

model.coint5 <- lm(pkn5$Close ~ lts5$Close)
summary(model.coint5)
## 
## Call:
## lm(formula = pkn5$Close ~ lts5$Close)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -10.0379  -3.0494  -0.2781   3.0094  11.1425 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -9.09111    1.33600  -6.805 2.08e-11 ***
## lts5$Close   0.68730    0.00996  69.007  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 4.098 on 745 degrees of freedom
## Multiple R-squared:  0.8647, Adjusted R-squared:  0.8645 
## F-statistic:  4762 on 1 and 745 DF,  p-value: < 2.2e-16

5. Testuję stacjonarność reszt

testdf(variable = model.coint5$residuals, adf_order = 1)

##   order       adf p_adf   bgodfrey        p_bg
## 2     0 -6.118491 <1pct 8.52785257 0.003497521
## 3     1 -5.272414 <1pct 0.03817377 0.845094628

W celu zbadania kointegracji obu zmiennych należy oszacować wektor kointegrujący i przetestować stacjonarność reszt.

Wartość statystyki DF równa -6.118491 oznacza brak podstaw do odrzucenia hipotezy zerowej o stacjonarnosci reszt.

Wartość statystyki BG równa 8.52785257 oraz wartość p-value na poziomie 0.003497521 oznaczają brak podstaw do odrzucenia hipotezy o braku autokorelacji reszt.

Oznacza to, że reszty z oszacowanego modelu są stacjonarne, co świadczy o kointegracji spółek ING i Pekao SA.

6. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Skoro ceny obu spółek są skointegrowane, będę mógł zastosować strategię pair-trading. W tym celu muszę wyznaczyć sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd

model.coint5$czy_sd<-0
model.coint5$czy_zero<-0
for(i in 1:747){
  model.coint5$czy_sd[i] = ifelse(abs(model.coint5$residuals[i])>sd(model.coint5$residuals),1,0)
}

Sprawdzam, czy dwie kolejne wartości przecięły wartość 0.

for(i in 1:746){
  model.coint5$czy_zero[i+1]=ifelse(model.coint5$residuals[i]*model.coint5$residuals[i+1]>0,0,1)
}

model.coint5$zamk<-0
model.coint5$otw<-0
i=1
k=1

for(k in 1:length(model.coint5$residuals)){
  model.coint5$zamk[k]<-NA
}
for(k in 1:length(model.coint5$residuals)){
  model.coint5$otw[k]<-NA
}

Uzupełniam kolumny otw oraz zamk sygnałami.

while(i <= 747){
  if(model.coint5$czy_sd[i]==1){
    model.coint5$otw[i]=ifelse(model.coint5$residuals[i]>sd(model.coint5$residuals),sd(model.coint5$residuals),-sd(model.coint5$residuals))
    i<-i+1
    while(i <= 747){
      if(model.coint5$czy_zero[i]==1){
        model.coint5$zamk[i]=0
        break
      }
      i<-i+1
    }
  }
  i<-i+1
}

7. Wykres pozycji dla obu spółek

plot(pkn5$Date, model.coint5$residuals, type="l")
abline(h=c(sd(model.coint5$residuals),-sd(model.coint5$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pkn5$Date, y=model.coint5$otw, col="red", pch=17)
points(x=pkn5$Date, y=model.coint5$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (in-sample)", cex.main=0.85)

Powyższy wykres pokazuje sygnały otwarcia/zamknięcia pozycji w strategii pair-traiding w okresie in-sample. Opierają się one na wykresie wartości nierównowagi krótkookresowej, czyli tak naprawdę na wartościach reszt z oszacowanego modelu. Moment wychylenia się nierównowagi krótkookresowej poza wartość jednego odchylenia standardowego (czerwona linia) stanowi sygnał otwarcia pozycji (czerwony trójkąt). Pozycję zamykam (niebieski kwadrat), kiedy wartość reszty wróci do średniej, czyli wartości 0 (niebieska linia).

8. Obliczam zwroty dla obu spółek.

Zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

i=1
for(i in 2:747){
  model.coint5$zwrot_A[i] = pkn5$Open[i]/pkn5$Open[i-1]
}

i=1
for(i in 2:747){
  model.coint5$zwrot_B[i] = lts5$Open[i]/lts5$Open[i-1]
}

model.coint5$otw[model.coint5$otw==sd(model.coint5$residuals)]<-1
model.coint5$otw[model.coint5$otw==-sd(model.coint5$residuals)]<-1
model.coint5$otw[is.na(model.coint5$otw)]<-0
model.coint5$zamk[model.coint5$zamk==0]<-1
model.coint5$zamk[is.na(model.coint5$zamk)]<-0

9. Tworzę kolumnę pozycja.

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1
while(i<=746){
  repeat{
    model.coint5$pozycja[i]=0
    i<-i+1
    if(model.coint5$otw[i-1]!=0||i>746){break}
  }
  while(model.coint5$zamk[i]!=1&&i<=746){
    model.coint5$pozycja[i]=ifelse(model.coint5$residuals[i]>0,1,-1)
    i<-i+1
  }
  model.coint5$pozycja[i]=ifelse(model.coint5$residuals[i-1]>0,1,-1)
  i<-i+1
}

10. Obliczam zwrot z pozycji.

Kolumna zwrot z pozycji mówi, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
model.coint5$zwrot_z_pozycji[1]=0
for(i in 2:747){
  model.coint5$zwrot_z_pozycji[i] = model.coint5$pozycja[i]*(model.coint5$zwrot_B[i]-model.coint5$zwrot_A[i])
}

i=2
for(i in 2:747){
  if(model.coint5$otw[i]==1){
    model.coint5$zwrot_z_pozycji[i+1] = model.coint5$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:747){
  if(model.coint5$zamk[i]==1){
    model.coint5$zwrot_z_pozycji[i] = model.coint5$zwrot_z_pozycji[i]*0.999
  }
}

11. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

model.coint5$portfel[1]=1
i=2
for(i in 2:747){
  model.coint5$portfel[i] = model.coint5$zwrot_z_pozycji[i-1]*model.coint5$portfel[i-1]+model.coint5$portfel[i-1]
}

12. Obliczam dzienne zwroty z indeksu WIG20

Liczę też zwrot z portfela dla strategii buy&hold dla WIG20.

i=2
wig20$zwrot<-0
for(i in 2:747){
  wig20$zwrot[i] = wig20$Open[i]/wig20$Open[i-1]
}

wig20$portfel[1]=1
i=2
for(i in 2:747){
  wig20$portfel[i] = wig20$zwrot[i]*wig20$portfel[i-1]
}

13. Wykres equity line

plot(pkn5$Date, model.coint5$portfel, type="l", ylim=c(0.5, 1.5), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pkn5$Date, wig20$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "topleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (in-sample)", cex.main=.85)

Praktycznie do połowy roku 2011 widoczna jest wyraźna przewaga zyskowności strategii opartej na kointegracji nad buy&hold na WIG20. Jej rentowność stabilnie utrzymuje się ponad rentownością strategii buy&hold na indeks WIG20.

Okres out-of-sample

1. DANE

Importuję dane i dokonuję niezbędnych poprawek na tych samych spółkach w okresie out-of-sample.

url1 <- "http://stooq.com/q/d/l/?s=ing&i=d"
pknOOS5 <- read.csv(url1,
                   header = TRUE,
                   sep = ",",
                   dec = ".",
                   stringsAsFactors = F)
pknOOS5$Date <- as.Date(pknOOS5$Date)
pknOOS5 <- pknOOS5[, c("Date","Open", "Close")]

url2 <- "http://stooq.com/q/d/l/?s=peo&i=d"
ltsOOS5 <- read.csv(url2,
                   header = TRUE,
                   sep = ",",
                   dec = ".",
                   stringsAsFactors = F)
ltsOOS5$Date <- as.Date(ltsOOS5$Date)
ltsOOS5 <- ltsOOS5[, c("Date", "Open", "Close")]

url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20OOS <- read.csv(url1,
                     header = TRUE,
                     sep = ",",
                     dec = ".",
                     stringsAsFactors = F)
wig20OOS$Date <- as.Date(wig20OOS$Date)
wig20OOS <- wig20OOS[, c("Date","Open", "Close")]

pknOOS5<-pknOOS5[!(pknOOS5$Date < "2014-01-01"), ]
pknOOS5<-pknOOS5[!(pknOOS5$Date > "2015-06-30"), ]

ltsOOS5<-ltsOOS5[!(ltsOOS5$Date < "2014-01-01"), ]
ltsOOS5<-ltsOOS5[!(ltsOOS5$Date > "2015-06-30"), ]

wig20OOS<-wig20OOS[!(wig20OOS$Date < "2014-01-01"), ]
wig20OOS<-wig20OOS[!(wig20OOS$Date > "2015-06-30"), ]

pknOOS5$Date <- as.Date(pknOOS5$Date, format = "%Y-%m-%d")
ltsOOS5$Date <- as.Date(ltsOOS5$Date, format = "%Y-%m-%d")
wig20OOS$Date <- as.Date(wig20OOS$Date, format = "%Y-%m-%d")

2. Wykres cen zamknięcia:

plot(pknOOS5$Date, pknOOS5$Close, type="l", ylim=c(100,190), ylab = "Cena zamknięcia", xlab = "Data", col="lightblue3")
lines(pknOOS5$Date, ltsOOS5$Close, col="firebrick3")
title("Ceny zamknięcia ING i Pekao (out-of-sample)", cex.main=0.85)

3. Obliczam reszty.

Obliczeń dokonuję na podstawie wektora kointegrującego oszacowanego dla okresu in-sample.

pknOOS5$teor=model.coint5$coefficients[1]+model.coint5$coefficients[2]*ltsOOS5$Close
pknOOS5$reszty=pknOOS5$Close-pknOOS5$teor

4. Wyznaczam sygnały otwarcia oraz zamknięcia pozycji.

Sprawdzam, czy dana reszta przekroczyła wartość jednego odchylenia standardowego i zaznaczam to w kolumnie czy_sd. Sprawdzam też, czy dwie kolejne wartości przecięły wartość 0.

pknOOS5$czy_sd<-0
pknOOS5$czy_zero<-0

for(i in 1:372){
  pknOOS5$czy_sd[i] = ifelse(abs(pknOOS5$reszty[i])>sd(model.coint5$residuals),1,0)
}

for(i in 1:371){
  pknOOS5$czy_zero[i+1]=ifelse(pknOOS5$reszty[i]*pknOOS5$reszty[i+1]>0,0,1)
}

pknOOS5$zamk<-0
pknOOS5$otw<-0
i=1
k=1

for(k in 1:length(pknOOS5$reszty)){
  pknOOS5$zamk[k]<-NA
}
for(k in 1:length(pknOOS5$reszty)){
  pknOOS5$otw[k]<-NA
}

while(i <= 372){
  if(pknOOS5$czy_sd[i]==1){
    pknOOS5$otw[i]=ifelse(pknOOS5$reszty[i]>sd(model.coint5$residuals),sd(model.coint5$residuals),-sd(model.coint5$residuals))
    i<-i+1
    while(i <= 372){
      if(pknOOS5$czy_zero[i]==1){
        pknOOS5$zamk[i]=0
        break
      }
      i<-i+1
    }
  }
  i<-i+1
}

5. Wykres pozycji dla obu spółek

plot(pknOOS5$Date, pknOOS5$reszty, type="l")
abline(h=c(sd(model.coint5$residuals),-sd(model.coint5$residuals)), col="red")
abline(h=0, col = "blue")
points(x=pknOOS5$Date, y=pknOOS5$otw, col="red", pch=17)
points(x=pknOOS5$Date, y=pknOOS5$zamk, col="blue", pch=15)
title("Sygnały otwarcia/zamknięcia (out-of-sample)", cex.main=0.85)

6. Obliczam zwroty dla obu spółek.

Następnie zmieniam wartości w kolumnach otw oraz zamk (służące do przygotowania wykresu pozycji) na wartości 1 oznaczające sygnały.

i=1
for(i in 2:372){
  pknOOS5$zwrot_A[i] = pknOOS5$Open[i]/pknOOS5$Open[i-1]
}

i=1
for(i in 2:372){
  pknOOS5$zwrot_B[i] = ltsOOS5$Open[i]/ltsOOS5$Open[i-1]
}

pknOOS5$otw[pknOOS5$otw==sd(model.coint5$residuals)]<-1
pknOOS5$otw[pknOOS5$otw==-sd(model.coint5$residuals)]<-1
pknOOS5$otw[is.na(pknOOS5$otw)]<-0
pknOOS5$zamk[pknOOS5$zamk==0]<-1
pknOOS5$zamk[is.na(pknOOS5$zamk)]<-0

7. Tworzę kolumnę pozycja

Przyjmuje ona 3 wartości:

  1. 0 - w przypadku braku otwartych pozycji
  2. 1 - w przypadku pozycji długiej w spółce B i krótkiej w spółce A
  3. -1 - w przypadku pozycji krótkiej w spółce B i długiej w spółce A
i=1
while(i<=371){
  repeat{
    pknOOS5$pozycja[i]=0
    i<-i+1
    if(pknOOS5$otw[i-1]!=0||i>371){break}
  }
  while(pknOOS5$zamk[i]!=1&&i<=371){
    pknOOS5$pozycja[i]=ifelse(pknOOS5$reszty[i]>0,1,-1)
    i<-i+1
  }
  pknOOS5$pozycja[i]=ifelse(pknOOS5$reszty[i-1]>0,1,-1)
  i<-i+1
}

8. Obliczam zwrot z pozycji.

Mówi on, ile wynosi łaczny zwrot z pozycji długiej oraz krótkiej. Dokonuję korekty zwrotu o koszty transakcyjne.

i=2
pknOOS5$zwrot_z_pozycji[1]=0
for(i in 2:372){
  pknOOS5$zwrot_z_pozycji[i] = pknOOS5$pozycja[i]*(pknOOS5$zwrot_B[i]-pknOOS5$zwrot_A[i])
}

i=2
for(i in 2:372){
  if(pknOOS5$otw[i]==1){
    pknOOS5$zwrot_z_pozycji[i+1] = pknOOS5$zwrot_z_pozycji[i+1]*0.999
  }
}

i=2
for(i in 2:372){
  if(pknOOS5$zamk[i]==1){
    pknOOS5$zwrot_z_pozycji[i] = pknOOS5$zwrot_z_pozycji[i]*0.999
  }
}

9. Obliczam zwrot z portfela.

Posłuży on do wyznaczenia equity line.

pknOOS5$portfel[1]=1
i=2
for(i in 2:372){
  pknOOS5$portfel[i] = pknOOS5$zwrot_z_pozycji[i-1]*pknOOS5$portfel[i-1]+pknOOS5$portfel[i-1]
}

10. Obliczam dzienne zwroty z indeksu WIG20.

Liczę również wartość portfela dla strategii buy&hold dla WIG20.

i=1
wig20OOS$zwrot<-0
for(i in 2:372){
  wig20OOS$zwrot[i] = wig20OOS$Open[i]/wig20OOS$Open[i-1]
}

wig20OOS$portfel[1]=1
i=2
for(i in 2:372){
  wig20OOS$portfel[i] = wig20OOS$zwrot[i]*wig20OOS$portfel[i-1]
}

11. Wykres equity line.

plot(pknOOS5$Date, pknOOS5$portfel, type="l", ylim=c(0.6, 1.2), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pknOOS5$Date, wig20OOS$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "topleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (out-of-sample)", cex.main=.85)

12. Statystyki opisowe dla okresu in-sample

model.coint5$zwrot_z_pozycji_na<-model.coint5$zwrot_z_pozycji
model.coint5$zwrot_z_pozycji_na[model.coint5$zwrot_z_pozycji_na==0]<-NA
annual.return<-Return.annualized(na.omit(model.coint5$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sd<-sd.annualized(na.omit(model.coint5$zwrot_z_pozycji_na), scale = 252)
beta<-lm(model.coint5$portfel~wig20$portfel)
drawdown<-maxdrawdown(model.coint5$portfel)
wsp.sharpe<-sharpe(model.coint5$portfel)
in.sample<-c(annual.return, annual.sd, wsp.sharpe, beta$coefficients[2], drawdown$maxdrawdown)

13. Statystyki opisowe dla okresu out-of-sample.

pknOOS5$zwrot_z_pozycji_na<-pknOOS5$zwrot_z_pozycji
pknOOS5$zwrot_z_pozycji_na[pknOOS5$zwrot_z_pozycji_na==0]<-NA
annual.returnOOS<-Return.annualized(na.omit(pknOOS5$zwrot_z_pozycji_na), scale = 252, geometric = FALSE)
annual.sdOOS<-sd.annualized(na.omit(pknOOS5$zwrot_z_pozycji_na), scale = 252)
betaOOS<-lm(pknOOS5$portfel~wig20OOS$portfel)
drawdownOOS<-maxdrawdown(pknOOS5$portfel)
wsp.sharpeOOS<-sharpe(pknOOS5$portfel)
wskaźnik<-c("zannualizowana stopa zwrotu", "zannualizowane odchylenie standardowe", "współczynnik Sharpe'a", "współczynnik beta", "maksymalne obsunięcie kapitału")
outof.sample<-c(annual.returnOOS, annual.sdOOS, wsp.sharpeOOS, betaOOS$coefficients[2], drawdownOOS$maxdrawdown)
statystyki<-data.frame(wskaźnik, in.sample, outof.sample)
knitr::kable(statystyki)
wskaźnik in.sample outof.sample
zannualizowana stopa zwrotu 0.1532631 -0.0781959
zannualizowane odchylenie standardowe 0.3154959 0.2966015
współczynnik Sharpe’a 0.2208379 -0.4603994
współczynnik beta -0.3650158 -0.1064613
maksymalne obsunięcie kapitału 0.1650844 0.3036168

Wstępna analiza wykresów equity line pokazuje, że otrzymane wyniki są zgodne z przewidywanymi. W okresie in-sample strategi oparta na kointegracji przynosi zyski w porównaniu z okresem out-of-sample.

Wartość zannualizowanej stopy zwrotu w okresie in-sample równa 0.1532631 mówi o tym, że strategia ta przyniosła zyski. Dla okresu out-of-sample ta wartość wynosi -0.0781959 co oznacza, że przyniosła straty.

Zannualizowane odchylenie standardowe dla okresu in-sample wyniosło 31.55%, dla out-of-sample 29.66%

Współczynnik Sharpe’a dla okresu in-sample wyniósł 0.2208379. Oznacza to, że zyskowność strategii była większa niż stopa wolna od ryzyka za ten okres. Przeciwne wnioski można wysnuć z wartości współczynnika Sharpe’a dla okresu out-of-sample równego -0.4603994.

Współczynnik beta wyniósł -0.3650158 w okresie in-sample, co świadczy o ujemnej korelacji strategii opartej na kointegracji względem buy-and-hold dla indeksu WIG20. Wartośc -0.1064613 dla okresu out-of-sample również świadczy o ujemnej korelacji względem indeksu WIG20.

Wskaźnik maksymalnego obsunięcia kapitału równy 0.1650844 w okresie in sample, w porównaniu do 0.3036168 świadczy o mniejszym ryzyku związanym ze strategią pair trading w okresie in-sample.

Portfel syntetyczny

Okres in-sample

1. Pobranie i wstępna obróbka notowań WIG20

wig20<-0
url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20 <- read.csv(url1,
                  header = TRUE,
                  sep = ",",
                  dec = ".",
                  stringsAsFactors = F)
wig20$Date <- as.Date(wig20$Date)
wig20 <- wig20[, c("Date","Open", "Close")]
wig20<-wig20[!(wig20$Date < "2011-01-01"), ]
wig20<-wig20[!(wig20$Date > "2013-12-31"), ]

#Przekształcam datę na właściwy format
wig20$Date <- as.Date(wig20$Date, format = "%Y-%m-%d")

2. Obliczenie dziennych zwrotów WIG20.

Obliczyłem też całkowity zwrot z portfela buy&hold.

i=1
wig20$zwrot<-0
for(i in 2:747){
  wig20$zwrot[i] = wig20$Open[i]/wig20$Open[i-1]
}

wig20$portfel[1]=1
i=2
for(i in 2:747){
  wig20$portfel[i] = wig20$zwrot[i]*wig20$portfel[i-1]
}

syntetyczny_in<-0
syntetyczny_in$portfel<-0
i=1

3. Obliczenie zwrotów jako średniej ze wszystkich analizowanych par spółek.

for(i in 1:747){
  syntetyczny_in$portfel[i]<-(model.coint$portfel[i]+model.coint2$portfel[i]+model.coint3$portfel[i]+model.coint4$portfel[i]+model.coint5$portfel[i])/5
}

syntetyczny_in$zwrot_z_pozycji<-0
i=1
for(i in 1:747){
  syntetyczny_in$zwrot_z_pozycji[i]<-(model.coint$zwrot_z_pozycji[i]+model.coint2$zwrot_z_pozycji[i]+model.coint3$zwrot_z_pozycji[i]+model.coint4$zwrot_z_pozycji[i]+model.coint5$zwrot_z_pozycji[i])/5
}

Okres out-of-sample

1. Import danych

wig20OOS<-0
url1 <- "http://stooq.com/q/d/l/?s=wig20&i=d"
wig20OOS <- read.csv(url1,
                     header = TRUE,
                     sep = ",",
                     dec = ".",
                     stringsAsFactors = F)
wig20OOS$Date <- as.Date(wig20OOS$Date)
wig20OOS<-wig20OOS[, c("Date","Open", "Close")]
wig20OOS<-wig20OOS[!(wig20OOS$Date < "2014-01-01"), ]
wig20OOS<-wig20OOS[!(wig20OOS$Date > "2015-06-30"), ]
wig20OOS$Date <- as.Date(wig20OOS$Date, format = "%Y-%m-%d")

2. Obliczenie dziennych zwrotów WIG20

i=1
wig20OOS$zwrot<-0
for(i in 2:372){
  wig20OOS$zwrot[i] = wig20OOS$Open[i]/wig20OOS$Open[i-1]
}

wig20OOS$portfel[1]=1
i=2
for(i in 2:372){
  wig20OOS$portfel[i] = wig20OOS$zwrot[i]*wig20OOS$portfel[i-1]
}

syntetyczny_out<-0
syntetyczny_out$portfel<-0
i=1

3. Obliczenie zwrotów jako średniej ze wszystkich analizowanych spółek

for(i in 1:372){
  syntetyczny_out$portfel[i]<-(pknOOS$portfel[i]+pknOOS2$portfel[i]+pknOOS3$portfel[i]+pknOOS4$portfel[i]+pknOOS5$portfel[i])/5
}

syntetyczny_out$zwrot_z_pozycji<-0
i=1
for(i in 1:372){
  syntetyczny_out$zwrot_z_pozycji[i]<-(pknOOS$zwrot_z_pozycji[i]+pknOOS2$zwrot_z_pozycji[i]+pknOOS3$zwrot_z_pozycji[i]+pknOOS4$zwrot_z_pozycji[i]+pknOOS5$zwrot_z_pozycji[i])/5
}

Wykresy wartości portfela i statystyki opisowe

plot(pkn2$Date, syntetyczny_in$portfel, type="l", ylim=c(0.5,2), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pkn2$Date, wig20$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "topleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (in-sample)", cex.main=.85)

plot(pknOOS2$Date, syntetyczny_out$portfel, type="l", ylim=c(0.8,1.1), ylab = "Zwrot", xlab = "Data", col="darkgoldenrod")
lines(pknOOS2$Date, wig20OOS$portfel, type="l", col="grey")
legend(bty="n", ncol=3, "topleft", "groups", c("Equity line","buy&hold"), lty=1, col=c("darkgoldenrod", "grey"), pt.cex=1.7, cex=0.6)
title("Equity line (out-of-sample)", cex.main=.85)

syntetyczny_in$zwrot_z_pozycji_na<-syntetyczny_in$zwrot_z_pozycji
syntetyczny_in$zwrot_z_pozycji_na[syntetyczny_in$zwrot_z_pozycji_na==0]<-NA
annual.return<-Return.annualized(na.omit(syntetyczny_in$zwrot_z_pozycji_na), scale = 252, geometric = TRUE)
annual.sd<-sd.annualized(na.omit(syntetyczny_in$zwrot_z_pozycji_na), scale = 252)
beta<-lm(syntetyczny_in$portfel~wig20$portfel)
drawdown<-maxdrawdown(syntetyczny_in$portfel)
wsp.sharpe<-sharpe(syntetyczny_in$portfel)
in.sample<-c(annual.return, annual.sd, wsp.sharpe, beta$coefficients[2], drawdown$maxdrawdown)
syntetyczny_out$zwrot_z_pozycji_na<-syntetyczny_out$zwrot_z_pozycji
syntetyczny_out$zwrot_z_pozycji_na[syntetyczny_out$zwrot_z_pozycji_na==0]<-NA
annual.returnOOS<-Return.annualized(na.omit(syntetyczny_out$zwrot_z_pozycji_na), scale = 252, geometric = TRUE)
annual.sdOOS<-sd.annualized(na.omit(syntetyczny_out$zwrot_z_pozycji_na), scale = 252)
betaOOS<-lm(syntetyczny_out$portfel~wig20OOS$portfel)
drawdownOOS<-maxdrawdown(syntetyczny_out$portfel)
wsp.sharpeOOS<-sharpe(syntetyczny_out$portfel)
wskaźnik<-c("zannualizowana stopa zwrotu", "zannualizowane odchylenie standardowe", "współczynnik Sharpe'a", "współczynnik beta", "maksymalne obsunięcie kapitału")
outof.sample<-c(annual.returnOOS, annual.sdOOS, wsp.sharpeOOS, betaOOS$coefficients[2], drawdownOOS$maxdrawdown)
statystyki<-data.frame(wskaźnik, in.sample, outof.sample)
knitr::kable(statystyki)
wskaźnik in.sample outof.sample
zannualizowana stopa zwrotu 0.1253083 -0.0494932
zannualizowane odchylenie standardowe 0.1179321 0.1445631
współczynnik Sharpe’a 0.8707299 -0.3533116
współczynnik beta -1.0141314 0.6974272
maksymalne obsunięcie kapitału 0.1651575 0.2190049

Analiza statystyk opisowych portfela syntetycznego zbudowanego na podstawie wszystkich 5 analizowanych par spółek potwierdza, że strategia pair trading oparta na kointegracji przynosi zyski w okresie in-sample, lecz wiąże się ze stratami w okresie out-of-sample. Można to również odczytać z wykresów equity curve.

Zannualizowana stopa zwrotu na poziomie 0.1253083 dla okresu in-sample wobec -0.0494932 w okresie out-of sample potwierdza tę obserwację.

W obu wykresach zannualizowane odchylenie standardowe pozostawało w granicach 11-14%.

Współczynnik Sharpe’a wyniósł 0.8707299 w okresie in-sample co oznacza, że portfel syntetyczny był w tym okresie bardziej zyskowny niż stopa wolna od ryzyka. Przeciwny wniosek wiąże się z wartością -0.3533116 za okres out-of-sample.

Współczynnik beta dla okresu in-sample na poziomie -1.0141314 pokazuje silną ujemną korelację portfela względem strategii buy&hold dla WIG20. W okresie out-of-sample beta wyniosła 0.6974272 co z kolei świadczy o dodatniej korelacji.

Wskaźnik maksymalnego obsunięcia kapitału na poziomie 0.1651575 w okresie in-sample w porównaniu do 0.2190049 w out-of-sample świadczy o większym ryzyku związanym z portfelem w okresie out-of-sample

Podsumowanie

Okazuje się, że znalazły się takie skointegrowane pary spółek, które dawały przeciwne do spodziewanych rezultatów (lepsze wyniki w okresie out-of-sample niż in-sample). Wiąże się to często z liczbą sygnałów a w konsekwencji wykonanych transakcji. Bywa, że relacja między badanymi spółkami jest różna w obu okresach i nie można wtedy skutecznie wykorzystać tej strategii dla okresu out-of-sample. Skuteczność badanej strategii należałoby sprawdzić poprzez zbadanie większej liczby par spółek. Trudno jest wyciągnąć jednoznaczne wnioski dotyczące pair tradingu opartego o kointegrację na tak małej próbie, jednak wyniki portfela syntetycznego pokazują, że nawet jeżeli wynik jednej pary spółek okaże się być odwrotny od oczekiwanego, portfel zbudowany na kilku takich parach okazuje się przynosić zyski.