#Data sets are imported into R

library(readxl)
BEF_Thesis_Volkov_Fisher_Hedging_BTC_FX <- read_excel("C:/Users/volkoa3/Desktop/Thesis/BEF_Thesis_Volkov_Fisher&Hedging_BTC_FX.xlsx")
View(BEF_Thesis_Volkov_Fisher_Hedging_BTC_FX)
BEF_Thesis_Volkov_Inflation <- read_excel("C:/Users/volkoa3/Desktop/Thesis/BEF_Thesis_Volkov_Inflation.xlsx")
View(BEF_Thesis_Volkov_Inflation)
BEF_Thesis_Volkov_Gold_SP500 <- read_excel("C:/Users/volkoa3/Desktop/Thesis/BEF_Thesis_Volkov_Gold_SP500.xlsx")
View(BEF_Thesis_Volkov_Gold_SP500)
BEF_Thesis_Volkov_TIPS <- read_excel("C:/Users/volkoa3/Desktop/Thesis/BEF_Thesis_Volkov_TIPS.xlsx")
View(BEF_Thesis_Volkov_TIPS)
BEF_Thesis_Volkov_VAR_BTC <- read_excel("C:/Users/volkoa3/Desktop/Thesis/BEF_Thesis_Volkov_VAR_BTC.xlsx")
View(BEF_Thesis_Volkov_VAR_BTC)


#Setting Bitcoin prices

BTCUSD=BEF_Thesis_Volkov_Fisher_Hedging_BTC_FX$`BTC USD`
BTCEUR=BTCUSD*BEF_Thesis_Volkov_Fisher_Hedging_BTC_FX$`FX EUR`
BTCCZK=BTCUSD*BEF_Thesis_Volkov_Fisher_Hedging_BTC_FX$`FX CZK`

#Setting Bitcoin returns

R_BTCUSD=diff(BTCUSD)/BTCUSD[-length(BTCUSD)]
R_BTCEUR=diff(BTCEUR)/BTCEUR[-length(BTCEUR)]
R_BTCCZK=diff(BTCCZK)/BTCCZK[-length(BTCCZK)]

#Setting inflation rates

USCPI=BEF_Thesis_Volkov_Inflation$`US CPI`
EUHICP=BEF_Thesis_Volkov_Inflation$`EU HICP`
CZHICP=BEF_Thesis_Volkov_Inflation$`CZ HICP`

#Setting gold price and returns

Gold=BEF_Thesis_Volkov_Gold_SP500$Gold
R_gold=diff(Gold)/Gold[-length(Gold)]

#Setting S&P500

SP500=BEF_Thesis_Volkov_Gold_SP500$SP500
R_SP500=diff(SP500)/SP500[-length(SP500)]

#Setting TIPS

TIPS=BEF_Thesis_Volkov_TIPS$TIPS

#Basic descriptive statistics for Bitcoin
library(psych)

min(R_BTCUSD)
min(R_BTCEUR)
min(R_BTCCZK)

max(R_BTCUSD)
max(R_BTCEUR)
max(R_BTCCZK)

mean(R_BTCUSD)
mean(R_BTCEUR)
mean(R_BTCCZK)

median(R_BTCUSD)
median(R_BTCEUR)
median(R_BTCCZK)

sd(R_BTCUSD)
sd(R_BTCEUR)
sd(R_BTCCZK)

skew(R_BTCUSD)
skew(R_BTCEUR)
skew(R_BTCCZK)

kurtosi(R_BTCUSD)
kurtosi(R_BTCEUR)
kurtosi(R_BTCCZK)

#Basic descriptive statistics for inflation

min(USCPI)
min(EUHICP)
min(CZHICP)

max(USCPI)
max(EUHICP)
max(CZHICP)

mean(USCPI)
mean(EUHICP)
mean(CZHICP)

median(USCPI)
median(EUHICP)
median(CZHICP)

sd(USCPI)
sd(EUHICP)
sd(CZHICP)

skew(USCPI)
skew(EUHICP)
skew(CZHICP)

kurtosi(USCPI)
kurtosi(EUHICP)
kurtosi(CZHICP)

#Basic descriptive statistics for gold, S&P500 and TIPS

min(R_gold)
min(R_SP500)
min(TIPS)

max(R_gold)
max(R_SP500)
max(TIPS)

mean(R_gold)
mean(R_SP500)
mean(TIPS)

median(R_gold)
median(R_SP500)
median(TIPS)

sd(R_gold)
sd(R_SP500)
sd(TIPS)

skew(R_gold)
skew(R_SP500)
skew(TIPS)

kurtosi(R_gold)
kurtosi(R_SP500)
kurtosi(TIPS)

#OLS model BTCUSD

Model_BTCUSD=lm(R_BTCUSD ~ USCPI)
summary(Model_BTCUSD)

#OLS model BTCEUR

Model_BTCEUR=lm(R_BTCEUR ~ EUHICP)
summary(Model_BTCEUR)

#OLS model BTCCZK

Model_BTCCZK=lm(R_BTCCZK ~ CZHICP)
summary(Model_BTCCZK)

#Constructing CI 95%

confint(Model_BTCUSD,level = 0.95)
confint(Model_BTCEUR,level = 0.95)
confint(Model_BTCCZK,level = 0.95)

#Calculting hedging demand for US

corUS=cor(R_BTCUSD,USCPI)
varUSCPI=var(USCPI)
varBTCUSD=var(R_BTCUSD)
HDUS=corUS*sqrt(varUSCPI/varBTCUSD)

#Calculting hedging demand for EU

corEU=cor(R_BTCEUR,EUHICP)
varEUHICP=var(EUHICP)
varBTCEUR=var(R_BTCEUR)
HDEU=corEU*sqrt(varEUHICP/varBTCEUR)

#Calculting hedging demand for CZ

corCZ=cor(R_BTCCZK,CZHICP)
varCZHICP=var(CZHICP)
varBTCCZK=var(R_BTCCZK)
HDCZ=corCZ*sqrt(varCZHICP/varBTCCZK)

#Building VAR Model

library(vars)

tsBTC=ts(R_BTCUSD,start = c(2014,11,frequency(12)),end=c(2022,10,frequency(12)),frequency = 12)
tsCPI=ts(USCPI,start = c(2014,11,frequency(12)),end=c(2022,10,frequency(12)),frequency = 12)
tsGOLD=ts(R_gold,start = c(2014,11,frequency(12)),end=c(2022,10,frequency(12)),frequency = 12)
tsSP500=ts(R_SP500,start = c(2014,11,frequency(12)),end=c(2022,10,frequency(12)),frequency = 12)
tsTIPS=ts(TIPS,start = c(2014,11,frequency(12)),end=c(2022,10,frequency(12)),frequency = 12)

varlags=cbind(tsBTC,tsCPI,tsGOLD,tsSP500,tsTIPS)
lagselect=VARselect(varlags,lag.max=3, type="const")
lagselect$selection

ModelVAR=VAR(varlags, p=1, type="const", season=NULL, exog=NULL)
summary(ModelVAR)

BTC_n=BEF_Thesis_Volkov_VAR_BTC$`BTC n`
BTC_n1=BEF_Thesis_Volkov_VAR_BTC$`BTC n-1`
CPI_n1=BEF_Thesis_Volkov_VAR_BTC$`CPI n-1`
GOLD_n1=BEF_Thesis_Volkov_VAR_BTC$`Gold n-1`
SP500_n1=BEF_Thesis_Volkov_VAR_BTC$`SP500 n-1`
TIPS_n1=BEF_Thesis_Volkov_VAR_BTC$`TIPS n-1`

Model_BTCn=lm( BTC_n ~ BTC_n1 + CPI_n1 + GOLD_n1 + SP500_n1 +TIPS_n1 )
summary(Model_BTCn)

#Heteroskedasticity and autocorrelation tests

bptest(Model_BTCUSD)
bptest(Model_BTCEUR)
bptest(Model_BTCCZK)
bptest(Model_BTCn)

dwtest(Model_BTCUSD)
dwtest(Model_BTCEUR)
dwtest(Model_BTCCZK)
dwtest(Model_BTCn)

#Plotting graphs

library(ggplot2)
library(scales)

date1=BEF_Thesis_Volkov_Fisher_Hedging_BTC_FX$Date
date2=BEF_Thesis_Volkov_Inflation$Date

dfBTCUSD=data.frame(x=as.Date(date1),y=BTCUSD)
ggpBTCUSD=ggplot(dfBTCUSD, aes(x,y)) + geom_line(size=1.3, col="dodgerblue2")

ggpBTCUSD + scale_x_date(limits=as.Date(c("2014-10-01", "2022-10-01")),labels = date_format('%Y'))+
  scale_y_continuous(limits = c(0,70000), breaks = c(0,10000,20000,30000,40000,50000,60000,70000),labels=dollar_format(suffix = " $",prefix=""))+
  xlab("Year") +
  ylab("Bitcoin price")+
  theme_minimal()+
  theme(axis.title.y = element_text(vjust=3)) +
  theme(axis.title.x = element_text(vjust=3))

dfBTCEUR=data.frame(x=as.Date(date1),y=BTCEUR)
ggpBTCEUR=ggplot(dfBTCEUR, aes(x,y)) + geom_line(size=1.3, col="dodgerblue2")

ggpBTCEUR + scale_x_date(limits=as.Date(c("2014-10-01", "2022-10-01")),labels = date_format('%Y'))+
  scale_y_continuous(limits = c(0,60000), breaks = c(0,10000,20000,30000,40000,50000,60000),labels = dollar_format(suffix = " €",prefix=""))+
  xlab("Year") +
  ylab("Bitcoin price")+
  theme_minimal()+
  theme(axis.title.y = element_text(vjust=3)) +
  theme(axis.title.x = element_text(vjust=3))

dfBTCCZK=data.frame(x=as.Date(date1),y=BTCCZK)
ggpBTCCZK=ggplot(dfBTCCZK, aes(x,y)) + geom_line(size=1.3, col="dodgerblue2")

ggpBTCCZK + scale_x_date(limits=as.Date(c("2014-10-01", "2022-10-01")),labels = date_format('%Y'))+
  scale_y_continuous(limits = c(0,1400000), breaks = c(0,200000,400000,600000,800000,1000000,1200000,1400000),labels = dollar_format(suffix = " Kč",prefix=""))+
  xlab("Year") +
  ylab("Bitcoin price")+
  theme_minimal()+
  theme(axis.title.y = element_text(vjust=3)) +
  theme(axis.title.x = element_text(vjust=3))

dfRBTCUSD=data.frame(x=R_BTCUSD)
ggpRBTCUSD=ggplot(dfRBTCUSD, aes(x)) + geom_histogram(binwidth=0.1,color="black")
ggpRBTCUSD + geom_vline(aes(xintercept = mean(x)),color="red",linetype="dashed", size=2) +
  scale_y_continuous(limits=c(0,20),breaks = c (0,5,10,15,20)) +
  ylab("Frequency")+
  xlab("Bitcoin returns in USD")+
  theme_minimal()

dfRBTCEUR=data.frame(x=R_BTCEUR)
ggpRBTCEUR=ggplot(dfRBTCEUR,aes(x)) + geom_histogram(binwidth=0.09,color="black")
ggpRBTCEUR +  geom_vline(aes(xintercept = mean(x)),color="red",linetype="dashed", size=2) +
  scale_y_continuous(limits=c(0,20),breaks = c (0,5,10,15,20)) +
  ylab("Frequency")+
  xlab("Bitcoin returns in EUR")+
  theme_minimal()

dfRBTCCZK=data.frame(x=R_BTCCZK)
ggpRBTCCZK=ggplot(dfRBTCCZK, aes(x)) + geom_histogram(binwidth=0.1,color="black")
ggpRBTCCZK +  geom_vline(aes(xintercept = mean(x)),color="red",linetype="dashed", size=2) +
  scale_y_continuous(limits=c(0,20),breaks = c (0,5,10,15,20)) +
  ylab("Frequency")+
  xlab("Bitcoin returns in CZK")+
  theme_minimal()

dfInflation=data.frame(x=as.Date(date2),y1=USCPI,y2=EUHICP,y3=CZHICP)
ggpInflation=ggplot(dfInflation,aes(x)) + geom_line(aes(y=y1, color="US"))+          
  geom_line(aes(y=y2, color="Euro Area"))+
  geom_line(aes(y=y3,color="CZ"))

ggpInflation + scale_x_date(limits=as.Date(c("2014-11-01", "2022-10-01")),labels = date_format('%Y'))+
  scale_y_continuous(limits=c(-0.01,0.19),breaks = c(-0.01,0,0.01,0.03,0.05,0.07,0.09,0.11,0.13,0.15,0.17,0.19),labels=scales::percent)+
  xlab("Year") +
  ylab("Inflation rate")+
  theme_minimal()+
  theme(axis.title.y = element_text(vjust=3)) +
  theme(axis.title.x = element_text(vjust=3))+
  scale_color_manual(values=c("US"="dodgerblue2","Euro Area"="firebrick1","CZ"="black"))+
  theme(legend.position = "bottom")+
  theme(legend.title =element_blank())

dfUSCPI=data.frame(x=USCPI)
ggpUSCPI=ggplot(dfUSCPI,aes(x)) + geom_histogram(binwidth=0.01,color="black")
ggpUSCPI + 
  geom_vline(aes(xintercept = mean(x)),color="red",linetype="dashed", size=2) +
  scale_y_continuous(limits = c(0,40),breaks = c(0,10,20,30,40))+
  scale_x_continuous(breaks = c(-0.02,0,0.02,0.04,0.06,0.08))+
  ylab("Frequency")+
  xlab("Inflation rate US")+
  theme_minimal()

dfEUHICP=data.frame(x=EUHICP)
ggpEUHICP=ggplot(dfEUHICP,aes(x)) + geom_histogram(binwidth=0.01,color="black")
ggpEUHICP + geom_vline(aes(xintercept = mean(x)),color="red",linetype="dashed", size=2) +
  scale_y_continuous(limits = c(0,40),breaks = c(0,10,20,30,40))+
  scale_x_continuous(breaks = c(-0.02,0,0.02,0.04,0.06,0.08,0.10,0.12))+
  ylab("Frequency")+
  xlab("Inflation rate Euro Area")+
  theme_minimal()

dfCZHICP=data.frame(x=CZHICP)
ggpCZHICP=ggplot(dfCZHICP,aes(x)) + geom_histogram(binwidth=0.02,color="black")
ggpCZHICP + 
  geom_vline(aes(xintercept = mean(x)),color="red",linetype="dashed", size=2) +
  scale_y_continuous(limits = c(0,50),breaks = c(0,10,20,30,40,50))+
  scale_x_continuous(breaks = c(0,0.02,0.04,0.06,0.08,0.10,0.12,0.14,0.16,0.18))+
  ylab("Frequency")+
  xlab("Inflation rate CZ")+
  theme_minimal()

dfGold=data.frame(x=R_gold)
ggpgold=ggplot(dfGold,aes(x)) + geom_histogram(binwidth=0.02,color="black")
ggpgold +
  geom_vline(aes(xintercept = mean(x)),color="red",linetype="dashed", size=2) +
  scale_y_continuous(limits = c(0,30),breaks = c(0,10,20,30))+
  ylab("Frequency")+
  xlab("Gold returns")+
  theme_minimal()

dfSP500=data.frame(x=R_SP500)
ggpSP500=ggplot(dfSP500,aes(x)) + geom_histogram(binwidth=0.06,color="black")
ggpSP500 +
  geom_vline(aes(xintercept = mean(x)),color="red",linetype="dashed", size=2) +
  scale_y_continuous(limits = c(0,50),breaks = c(0,10,20,30,40,50))+
  ylab("Frequency")+
  xlab("S&P 500")+
  theme_minimal()

dfTIPS=data.frame(x=TIPS)
ggpTIPS=ggplot(dfTIPS,aes(x)) + geom_histogram(binwidth=0.005,color="black")
ggpTIPS +
  geom_vline(aes(xintercept = mean(x)),color="red",linetype="dashed", size=2) +
  scale_y_continuous(limits = c(0,40),breaks = c(0,10,20,30,40))+
  ylab("Frequency")+
  xlab("TIPS")+
  theme_minimal()

dfVAR=data.frame(x=as.Date(date2),y1=R_BTCUSD,y2=USCPI,y3=R_gold,y4=R_SP500,y5=TIPS)

ggpVAR=ggplot(dfVAR,aes(x))+
  geom_line(aes(y=y1,color="Bitcoin returns"))+
  geom_line(aes(y=y2,color="Inflation rate"))+
  geom_line(aes(y=y3,color="Gold returns"))+
  geom_line(aes(y=y4,color="S&P 500"))+
  geom_line(aes(y=y5,color="TIPS yield"))

ggpVAR +  
  scale_x_date(limits=as.Date(c("2014-11-01", "2022-10-01")),labels = date_format('%Y'))+
  xlab("Year") +
  ylab("Rate")+
  theme_minimal()+
  theme(axis.title.y = element_text(vjust=3)) +
  theme(axis.title.x = element_text(vjust=3))+
  scale_color_manual(values=c("Bitcoin returns"="dodgerblue2","Inflation rate"="firebrick1","Gold returns"="darkgoldenrod3","S&P 500"="chartreuse4","TIPS yield"="black"))+
  theme(legend.position = "bottom")+
  theme(legend.title =element_blank())

#Exporting models for Appendix

sink(file="BTCUSD.txt")
summary(Model_BTCUSD)
sink(file = NULL)

sink(file="BTCEUR.txt")
summary(Model_BTCEUR)
sink(file = NULL)

sink(file="BTCCZK.txt")
summary(Model_BTCCZK)
sink(file = NULL)

sink(file="VAR.txt")
summary(ModelVAR)
sink(file = NULL)