---
title: "<span style='font-size:25px;'>Assignment on - Quantitative Student Investment Fund</span>"
output:
  html_document: 
    code_download: true
    highlight: espresso
date: "`r Sys.Date()`"
editor_options: 
  markdown: 
    wrap: 72
---

###### **Submitted to:** **Jon Spinney**
###### **Submitted by:** Md Mahmudul Hasan
<style type="text/css"> body, td {font-size: 14px;} code.r{font-size: 12px;} pre {font-size: 14px} </style>
<script src="https://code.jquery.com/jquery-3.7.1.slim.min.js"></script>
<script type="text/javascript">
  $(document).ready(function() {
  $('a').attr('target', '_blank');
  });
</script>
<style>
.boxed {
  background: #232023;
  color: white;
  border: 0px solid #535353;
  margin: 0px auto;
  width: auto;
  padding: 10px;
  border-radius: 0px;
}
</style>

#  {.tabset}


```{r setup, include=FALSE}
knitr::opts_chunk$set(
	fig.align = "center",
	message = FALSE,
	warning = FALSE
)
library(fpp2)
library(tidyquant)
library(tidyverse)
library(tidyr)
library(quantmod)
library(dplyr)
library(ggthemes)
library(curl)
library(readxl)
library(cowplot)
library(ggpubr)
library(plotly)
library(ggplot2)
library(reshape2)
library(gt)
library(gtExtras)
library(knitr)
library(ggfortify)
library(moments)
library(skimr)
library(psych)
library(kableExtra)
library(forecast)
library(tseries)
library(tibble)
```

## Problem 1
<br>
<div class="boxed">
<br>
- Calculate the cumulative total price return for each of the 60 stocks as well as the TSX Composite Index from the beginning of the period to the end.
- Calculate the daily standard deviation of returns for each of the 60 stocks as well as the TSX Composite Index (annualize this number using the √N rule).
- Create a scatterplot in R plotting the standard deviation vs. the return for the year to date period for each of the 61 time series, with standard deviation on the x-axis. 
<br>
Use a different point style for the TSX Composite Index to highlight the index vs. the individual stocks.
<br>
</div>
<br>

```{r message=FALSE, warning=FALSE}
tsx_members <- read.csv("a1_tsx60memb.csv")
tsx_prices <- read.csv("a1_tsx60membpx.csv")
tsx_composite <- read.csv("a1_sptsxpx.csv")%>% 
  rename(Date = Dates, Price=S.P.TSX.Composite.Index.Level)

# Reshape the data to a long format
tsx_prices_long <- tsx_prices %>%
  pivot_longer(cols = -Date, names_to = "Ticker", values_to = "Price")

tsx_prices_long$Date <- as.Date(tsx_prices_long$Date)
tsx_composite$Date <- as.Date(tsx_composite$Date) 

# cumulative returns
cumulative_returns <- tsx_prices_long %>%
  group_by(Ticker) %>%
  summarize(cumulative_return = (last(Price) - first(Price)) / first(Price) * 100)

cumulative_returns %>%
  kbl(caption = "Cumulative Returns of S&P/TSX 60 Stocks") %>%
  kable_styling(full_width = F, position = "center", bootstrap_options = c("striped", "hover")) %>%
  scroll_box(height = "300px")  

```
```{r message=FALSE, warning=FALSE}
# daily returns and annualized standard deviation
daily_sd <- tsx_prices_long %>%
  arrange(Date) %>%
  group_by(Ticker) %>%
  mutate(daily_return = (Price - lag(Price)) / lag(Price)*100) %>%
  summarize(annualized_sd = sd(daily_return, na.rm = TRUE) * sqrt(252))

# Combine cumulative returns and standard deviations
tsx_data <- cumulative_returns %>%
  left_join(daily_sd, by = "Ticker")

# Now add the TSX Composite Index
tsx_composite_return <- (last(tsx_composite$Price) - first(tsx_composite$Price)) / first(tsx_composite$Price) * 100

# standard deviation for the TSX Composite Index
tsx_composite_sd <- tsx_composite %>%
  arrange(Date) %>%
  mutate(daily_return = (Price - lag(Price)) / lag(Price)*100) %>%
  summarize(annualized_sd = sd(daily_return, na.rm = TRUE) * sqrt(252))

# Combine stock data with the TSX Composite Index data
tsx_index_data <- data.frame(
  Ticker = "TSX Composite",
  cumulative_return = tsx_composite_return,
  annualized_sd = tsx_composite_sd$annualized_sd)
tsx_final_data <- bind_rows(tsx_data, tsx_index_data)

tsx_final_data %>%
  kbl(caption = "Combined Data of S&P/TSX Stocks and Index") %>%
  kable_styling(full_width = F, position = "center", bootstrap_options = c("striped", "hover")) %>%
  scroll_box(height = "300px")  

ggplot(tsx_final_data, aes(x = annualized_sd, y = cumulative_return, label = Ticker)) +
  geom_point(aes(shape = ifelse(Ticker == "TSX Composite", "Index", "Stock"), 
                 color = ifelse(Ticker == "TSX Composite", "Index", "Stock")), size = 3) +
  scale_shape_manual(values = c("Index" = 16, "Stock" = 17)) +  
  scale_color_manual(values = c("Index" = "#4584b6", "Stock" = "#ffde57")) +  
  labs(title = "Standard Deviation vs. Return for TSX Stocks and Composite Index",
       x = "Annualized Standard Deviation",
       y = "Cumulative Return (%)",
       color = "Legend",
       shape = "Legend") +
  theme_minimal(base_size = 15) +  
  geom_text(aes(label = Ticker), 
            hjust = 0.5, 
            vjust = -1,  
            size = 3, 
            check_overlap = TRUE) +  
  theme(legend.position = "top",  
        legend.box = "horizontal")  
```
<br>

## Problem 2
<br>
<div class="boxed">
<br>
In your *a1_tsx60memb* data frame, create a new column that assigns securities to a group based on market capitalization. 
<br>
Create three groups - one for securities with less than \$20B in market cap, one for securities with a market cap between \$20B and \$50B, and one for securities with a market cap above $50B. You may name the groups anything, but “Small”, “Mid”, and “Large” seems reasonable.
<br>
Use the *aggregate()* function in base R (or an equivalent function from the CRAN package data.table) to calculate the total weight of securities in each of your three market cap groups. Plot it in a bar chart.\
<br>
Create another bar chart but use the GICS Sectors as your grouping variable.
<br>
</div>
<br>

```{r message=FALSE, warning=FALSE}
tsx_members <- tsx_members %>%
  mutate(market_cap_group = case_when(
    Market.Cap < 20000000000 ~ "Small", 
    Market.Cap >= 20000000000 & Market.Cap <= 50000000000 ~ "Mid",  
    Market.Cap > 50000000000 ~ "Large" 
  ))
# total weight for each market cap group
market_cap_weight <- aggregate(Weight ~ market_cap_group, data = tsx_members, sum)

# bar chart
p <- ggplot(market_cap_weight, aes(x = market_cap_group, y = Weight)) +
  geom_bar(stat = "identity", fill = "#4584b6") +
  labs(title = "Total Weight by Market Cap Group", x = "Market Cap Group", y = "Total Weight") +
  theme_minimal()+
  theme(axis.text.x = element_blank())
interactive_plot <- ggplotly(p)%>%
  layout(hoverlabel = list(bgcolor = "#ffde57", font = list(color = "black")))
interactive_plot

# total weight for each GICS Sector
sector_weight <- aggregate(Weight ~ GICS.Sector, data = tsx_members, sum)

# a bar chart for GICS sectors
p_weight <- ggplot(sector_weight, aes(x = GICS.Sector, y = Weight)) +
  geom_bar(stat = "identity", fill = "#ffde57") +
  labs(title = "Total Weight by GICS Sector", x = "GICS Sector", y = "Total Weight") +
  theme_minimal() +
  theme(axis.text.x = element_blank())
  interactive_plot_weight <- ggplotly(p_weight)%>%
  layout(hoverlabel = list(bgcolor = "#4584b6", font = list(color = "white")))
interactive_plot_weight
```
<br>
<br>



## Problem 3
<br>
<div class="boxed">
<br>
The beta of a security is a risk measure that gauges the degree of market sensitivity of that security (high beta stocks have higher market risk, and vice versa) and is generally the slope of a regression line for the returns of that security regressed on the returns on the benchmark index.
<br>
For the securities in this sample, calculate the market beta using the TSX Composite Index as the benchmark (you will need to run 60 linear regressions and extract the slope coefficient from each). Use either a loop or preferably an apply() function to accomplish this.
<br>
Plot a barplot of the betas of the 60 stocks in your sample over the sample period. Also plot the R-Squared of the 60 regressions in a separate barplot. How effective is a one-factor model such as this in describing the returns?
<br>
</div>
<br>
<br>
```{r message=FALSE, warning=FALSE}
# daily returns for stocks 
tsx_prices_long <- tsx_prices_long %>%
  arrange(Date) %>%
  group_by(Ticker) %>%
  mutate(daily_return = (Price - lag(Price)) / lag(Price) * 100) %>%
  filter(!is.na(daily_return))  

# daily returns for TSX Composite Index
tsx_composite <- tsx_composite %>%
  arrange(Date) %>%
  mutate(daily_return = (Price - lag(Price)) / lag(Price) * 100) %>%
  filter(!is.na(daily_return))  

# A function to compute beta and R-squared for each stock
calculate_beta <- function(stock_returns, benchmark_returns) {
  model <- lm(stock_returns ~ benchmark_returns)
  beta <- coef(model)[2]  
  r_squared <- summary(model)$r.squared  
  return(c(beta = beta, r_squared = r_squared))}

# Data frame for beta and R-squared 
beta_results <- tsx_prices_long %>%
  group_by(Ticker) %>%
  summarize(
    beta_r_squared = list(calculate_beta(daily_return, tsx_composite$daily_return))
  ) %>%
  unnest_wider(beta_r_squared)

# Rename the columns for clarity
colnames(beta_results) <- c("Ticker", "Beta", "R_Squared")

# Plot the Betas
p_beta <- ggplot(beta_results, aes(x = reorder(Ticker, Beta), y = Beta)) +
  geom_bar(stat = "identity", fill = "#4584b6") +
  labs(title = "Betas of the 60 Stocks", x = "Stock Ticker", y = "Beta") +
  theme_minimal() +
  theme(axis.text.x = element_blank())

interactive_plot_beta <- ggplotly(p_beta)%>%
  layout(hoverlabel = list(bgcolor = "#ffde57", font = list(color = "black")))
interactive_plot_beta

# Plot the R-squared values
p_r_squared <- ggplot(beta_results, aes(x = reorder(Ticker, R_Squared), y = R_Squared)) +
  geom_bar(stat = "identity", fill = "#ffde57") +
  labs(title = "R-Squared of the 60 Regressions", x = "Stock Ticker", y = "R-Squared") +
  theme_minimal() +
  theme(axis.text.x = element_blank())
interactive_r_squared <- ggplotly(p_r_squared)%>%
  layout(hoverlabel = list(bgcolor = "#4584b6", font = list(color = "white")))
interactive_r_squared
```
<br>

## Problem 4
<br>
<div class="boxed">
<br>
<br>
Calculate the covariance matrix of the 60 securities in the TSX 60 Index, and use this data to estimate the volatility of:
<br>
• The equally weighted portfolio of all 60 securities,
• The the TSX Weighted portfolio (use the weights provided in the a1_tsx60memb.csv file)
• A portfolio of all 60 securities where the weight is inversely proportional to the beta you estimated in Problem 3 according to the following weight function:
$$\omega_i=\dfrac{\dfrac{1}{\beta_i}}{\sum_{i=1}^{N}\dfrac{1}{\beta_i}}$$
• An equally weighted portfolio containing only securities in the Financials GICS sector.
• Equally weighted portfolios for each of your three market capitalization groups.
Create a data frame containing the volatility for each of these 7 portfolios and plot them together in a barchart.
<br>
<br>
</div>


<br>

```{r message=FALSE, warning=FALSE}
# a pivot table of daily returns
tsx_returns_wide <- tsx_prices_long %>%
  select(Date, Ticker, daily_return) %>%
  pivot_wider(names_from = Ticker, values_from = daily_return)

# covariance matrix for the stocks (excluding the Date)
cov_matrix <- cov(tsx_returns_wide[, -1], use = "pairwise.complete.obs")
```


```{r message=FALSE, warning=FALSE}
# Equally weighted portfolio
equal_weights <- rep(1/60, 60)
equal_portfolio_variance <- t(equal_weights) %*% cov_matrix %*% equal_weights
equal_portfolio_sd <- sqrt(equal_portfolio_variance)

```


```{r message=FALSE, warning=FALSE}
# Extract the weights from the members data
tsx_weights <- tsx_members$Weight / 100 
tsx_portfolio_variance <- t(tsx_weights) %*% cov_matrix %*% tsx_weights
tsx_portfolio_sd <- sqrt(tsx_portfolio_variance)  
```


```{r message=FALSE, warning=FALSE}
# weights inversely proportional to beta
inverse_beta_weights <- 1 / beta_results$Beta
inverse_beta_weights <- inverse_beta_weights / sum(inverse_beta_weights)

# Inversely proportional to beta portfolio
inv_beta_portfolio_variance <- t(inverse_beta_weights) %*% cov_matrix %*% inverse_beta_weights
inv_beta_portfolio_sd <- sqrt(inv_beta_portfolio_variance)  
```


```{r message=FALSE, warning=FALSE}
# Filter for Financials sector
financials_tickers <- tsx_members %>%
  filter(GICS.Sector == "Financials") %>%
  pull(Ticker)
financials_tickers <- na.omit(financials_tickers)

# Subset the covariance matrix for Financials tickers
financials_cov_matrix <- cov_matrix[financials_tickers, financials_tickers]

# Equal weights for Financials sector (after removing NA)
financials_weights <- rep(1/length(financials_tickers), length(financials_tickers))

# Financials portfolio variance
financials_portfolio_variance <- t(financials_weights) %*% financials_cov_matrix %*% financials_weights
financials_portfolio_sd <- sqrt(financials_portfolio_variance)
```


```{r message=FALSE, warning=FALSE}
# Small Cap Portfolio
small_cap_tickers <- tsx_members %>%
  filter(market_cap_group == "Small") %>%
  pull(Ticker)

# Remove the invalid ticker from small_cap_tickers
small_cap_tickers <- small_cap_tickers[small_cap_tickers %in% rownames(cov_matrix)]

small_cap_weights <- rep(1/length(small_cap_tickers), length(small_cap_tickers))
small_cap_cov_matrix <- cov_matrix[small_cap_tickers, small_cap_tickers]
small_cap_portfolio_variance <- t(small_cap_weights) %*% small_cap_cov_matrix %*% small_cap_weights
small_cap_portfolio_sd <- sqrt(small_cap_portfolio_variance)

# Mid Cap Portfolio
mid_cap_tickers <- tsx_members %>%
  filter(market_cap_group == "Mid") %>%
  pull(Ticker)
mid_cap_weights <- rep(1/length(mid_cap_tickers), length(mid_cap_tickers))

# Remove the invalid tickers from mid_cap_tickers
mid_cap_tickers <- mid_cap_tickers[mid_cap_tickers %in% rownames(cov_matrix)]
mid_cap_cov_matrix <- cov_matrix[mid_cap_tickers, mid_cap_tickers]
mid_cap_weights <- rep(1 / length(mid_cap_tickers), length(mid_cap_tickers))
mid_cap_portfolio_variance <- t(mid_cap_weights) %*% mid_cap_cov_matrix %*% mid_cap_weights
mid_cap_portfolio_sd <- sqrt(mid_cap_portfolio_variance)

# Large Cap Portfolio
large_cap_tickers <- tsx_members %>%
  filter(market_cap_group == "Large") %>%
  pull(Ticker)
large_cap_weights <- rep(1/length(large_cap_tickers), length(large_cap_tickers))
large_cap_cov_matrix <- cov_matrix[large_cap_tickers, large_cap_tickers]
large_cap_portfolio_variance <- t(large_cap_weights) %*% large_cap_cov_matrix %*% large_cap_weights
large_cap_portfolio_sd <- sqrt(large_cap_portfolio_variance)
```


```{r message=FALSE, warning=FALSE}
# a data frame with the volatilities
volatility_data <- data.frame(
  Portfolio = c("Equally Weighted", "TSX Weighted", "Inverse Beta Weighted", 
                "Financials Sector", "Small Cap", "Mid Cap", "Large Cap"),
  Volatility = c(equal_portfolio_sd, tsx_portfolio_sd, inv_beta_portfolio_sd, 
                 financials_portfolio_sd, small_cap_portfolio_sd, mid_cap_portfolio_sd, large_cap_portfolio_sd))
```


```{r message=FALSE, warning=FALSE}
# ggplot bar chart
p <- ggplot(volatility_data, aes(x = Portfolio, y = Volatility)) +
  geom_bar(stat = "identity", fill = "#4584b6") +
  labs(title = "Volatility of Different Portfolios", x = "Portfolio", y = "Volatility (Standard Deviation)") +
  theme_minimal() +
  theme(axis.text.x = element_blank())
interactive_plot <- ggplotly(p) %>%
  layout(hoverlabel = list(bgcolor = "#ffde57", font = list(color = "black")))  
interactive_plot
```


## Problem 5
<br>
<div class="boxed">
<br>
<br>
Install the package quantmod and download prices of the S&P 500 using the function *getSymbols()* (ticker“ˆGSPC”). Calculate the daily returns and compare against a normal distribution using a *qqnorm()* and *qqline()* plot.
<br>
<br>
</div>
<br>
<br>
```{r message=FALSE, warning=FALSE}
# Download prices of the S&P 500 using getSymbols()
getSymbols("^GSPC", src = "yahoo", from = "2010-01-01", to = Sys.Date())

# daily returns- Using Adjusted Close prices to account for dividends and stock splits
sp500_prices <- Ad(GSPC)  
daily_returns <- diff(log(sp500_prices))  
daily_returns <- na.omit(daily_returns)  

# a Q-Q plot to compare against a normal distribution
qqnorm(daily_returns, main = "Q-Q Plot of Daily Returns of S&P 500")
qqline(daily_returns, col = "red")  
```

<br>

<br>


<br>

<span>  <span style="color: Steelblue;">**Md Mahmudul Hasan**</span> \
MQIM 3760573\
Faculty of Management \
University of New Brunswick \
mahmudul.hasan\@unb.ca \

<br>