## Import necessary library
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.4 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.1
## ✔ ggplot2 3.4.4 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.0
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(tidyquant)
## Loading required package: PerformanceAnalytics
## Loading required package: xts
## Loading required package: zoo
##
## Attaching package: 'zoo'
##
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
##
##
## ######################### Warning from 'xts' package ##########################
## # #
## # The dplyr lag() function breaks how base R's lag() function is supposed to #
## # work, which breaks lag(my_xts). Calls to lag(my_xts) that you type or #
## # source() into this session won't work correctly. #
## # #
## # Use stats::lag() to make sure you're not using dplyr::lag(), or you can add #
## # conflictRules('dplyr', exclude = 'lag') to your .Rprofile to stop #
## # dplyr from breaking base R's lag() function. #
## # #
## # Code in packages is not affected. It's protected by R's namespace mechanism #
## # Set `options(xts.warn_dplyr_breaks_lag = FALSE)` to suppress this warning. #
## # #
## ###############################################################################
##
## Attaching package: 'xts'
##
## The following objects are masked from 'package:dplyr':
##
## first, last
##
##
## Attaching package: 'PerformanceAnalytics'
##
## The following object is masked from 'package:graphics':
##
## legend
##
## Loading required package: quantmod
## Warning: package 'quantmod' was built under R version 4.3.3
## Loading required package: TTR
## Registered S3 method overwritten by 'quantmod':
## method from
## as.zoo.data.frame zoo
library(quantmod)
library(plotly)
##
## Attaching package: 'plotly'
##
## The following object is masked from 'package:ggplot2':
##
## last_plot
##
## The following object is masked from 'package:stats':
##
## filter
##
## The following object is masked from 'package:graphics':
##
## layout
library(DT)
## Warning: package 'DT' was built under R version 4.3.3
library(Metrics)
## Warning: package 'Metrics' was built under R version 4.3.3
library(lmtest)
## Warning: package 'lmtest' was built under R version 4.3.3
library(kableExtra)
##
## Attaching package: 'kableExtra'
##
## The following object is masked from 'package:dplyr':
##
## group_rows
library(quantmod)
# Download data
symbols <- c("MSFT", "AAPL", "INTC", "GOOG", "AMZN", "META", "UNH", "MCD", "ETN", "JPM", "V", "WMT", "XOM","MA", "PG", "ORCL", "COST", "JNJ", "HD", "MRK", "NKE", "ABBV", "NFLX", "CVX", "KO", "AMD", "QCOM", "ADBE", "PEP", "CRM")
getSymbols(symbols, from = "2019-05-31", to = "2024-05-18")
## [1] "MSFT" "AAPL" "INTC" "GOOG" "AMZN" "META" "UNH" "MCD" "ETN" "JPM"
## [11] "V" "WMT" "XOM" "MA" "PG" "ORCL" "COST" "JNJ" "HD" "MRK"
## [21] "NKE" "ABBV" "NFLX" "CVX" "KO" "AMD" "QCOM" "ADBE" "PEP" "CRM"
getSymbols("^GSPC", from = "2019-05-31", to = "2024-05-18") # S&P 500
## [1] "GSPC"
rf <- 0.0424 / 12 # Monthly risk free rate
# Calculate monthly log-return of each stock
monthly_returns <- function(symbol) {
monthly <- periodReturn(get(symbol), period = "monthly", type = "log")
colnames(monthly) <- symbol
return(monthly)
}
# Calculate monthly log-return of market
returns_data <- do.call(merge, lapply(symbols, monthly_returns))
market_returns <- periodReturn(GSPC, period = "monthly", type = "log")
colnames(market_returns) <- "Market"
# Calculate excess return of each stock
excess_returns <- returns_data - rf
# Calculate excess market return
market_excess_returns <- market_returns - rf
length(market_excess_returns)
## [1] 61
# Construct the Dataframe to store value of intercept coefficient and beta coefficient
betas <- list()
p_values_betas <- list()
alphas <- list()
p_value_alphas <- list()
# Using CAMP model to estimate beta coefficient of each stock (using linear regression)
for (symbol in symbols) {
model <- lm(excess_returns[, symbol] ~ market_excess_returns)
# Store alpha
alpha <- summary(model)$coefficients[1]
alphas[[symbol]] <- alpha
#Store P-Value of alpha
p_value_a <- summary(model)$coefficients["(Intercept)", "Pr(>|t|)"]
p_value_alphas[[symbol]] <- round(p_value_a, 4)
# Store beta
beta <- coef(model)["market_excess_returns"]
betas[[symbol]] <- beta
# Store P-Value of beta
p_value <- summary(model)$coefficients["market_excess_returns", "Pr(>|t|)"]
p_values_betas[[symbol]] <- round(p_value, 4)
}
# Combine all of data to a Dataframe
betas_df <- data.frame(
Symbol = symbols,
Alpha = unlist(alphas),
Beta = unlist(betas),
PValue_beta = unlist(p_values_betas),
PValue_alpha = unlist(p_value_alphas)
)
# Print Data Frame
print(betas_df)
## Symbol Alpha Beta PValue_beta PValue_alpha
## MSFT MSFT 0.0101192496 0.8862446 0.0000 0.0593
## AAPL AAPL 0.0116759996 1.2245193 0.0000 0.0947
## INTC INTC -0.0166839405 1.0895639 0.0000 0.1566
## GOOG GOOG 0.0082154804 1.0396709 0.0000 0.2541
## AMZN AMZN 0.0003027082 1.1266987 0.0000 0.9734
## META META 0.0037176481 1.1944103 0.0000 0.7839
## UNH UNH 0.0050677177 0.5726595 0.0000 0.4600
## MCD MCD -0.0033893929 0.7192985 0.0000 0.5155
## ETN ETN 0.0135260809 1.0186252 0.0000 0.0175
## JPM JPM -0.0007657541 1.1344930 0.0000 0.9094
## V V -0.0012712918 0.9563165 0.0000 0.8281
## WMT WMT 0.0035798416 0.4941423 0.0000 0.5466
## XOM XOM -0.0016810638 0.9329538 0.0000 0.8785
## MA MA -0.0014104659 1.0864504 0.0000 0.8381
## PG PG 0.0011782718 0.4169185 0.0004 0.8402
## ORCL ORCL 0.0037151505 1.0051208 0.0000 0.6185
## COST COST 0.0107026790 0.7685684 0.0000 0.1008
## JNJ JNJ -0.0046501256 0.5284799 0.0000 0.3983
## HD HD -0.0008454496 0.9927173 0.0000 0.8973
## MRK MRK 0.0026218085 0.4007275 0.0049 0.7165
## NKE NKE -0.0083129083 1.0534178 0.0000 0.3326
## ABBV ABBV 0.0049295568 0.5753174 0.0007 0.5592
## NFLX NFLX -0.0035636554 1.3400821 0.0000 0.8089
## CVX CVX -0.0055413896 1.0936172 0.0000 0.5602
## KO KO -0.0037630125 0.6105680 0.0000 0.4988
## AMD AMD 0.0138085295 1.6726159 0.0000 0.3639
## QCOM QCOM 0.0051306639 1.2687196 0.0000 0.5814
## ADBE ADBE -0.0033152247 1.2965920 0.0000 0.7225
## PEP PEP -0.0015627830 0.5347994 0.0000 0.7350
## CRM CRM -0.0023180603 1.2686408 0.0000 0.8197
# Taking expected market return
market_expected_return <- mean(market_returns)
# Calculate expected return of each stock by using CAPM model
expected_returns <- rf + betas_df$Beta * (market_expected_return - rf)
betas_df$ExpectedReturn <- expected_returns
print(betas_df)
## Symbol Alpha Beta PValue_beta PValue_alpha ExpectedReturn
## MSFT MSFT 0.0101192496 0.8862446 0.0000 0.0593 0.009858125
## AAPL AAPL 0.0116759996 1.2245193 0.0000 0.0947 0.012272263
## INTC INTC -0.0166839405 1.0895639 0.0000 0.1566 0.011309138
## GOOG GOOG 0.0082154804 1.0396709 0.0000 0.2541 0.010953070
## AMZN AMZN 0.0003027082 1.1266987 0.0000 0.9734 0.011574154
## META META 0.0037176481 1.1944103 0.0000 0.7839 0.012057386
## UNH UNH 0.0050677177 0.5726595 0.0000 0.4600 0.007620187
## MCD MCD -0.0033893929 0.7192985 0.0000 0.5155 0.008666694
## ETN ETN 0.0135260809 1.0186252 0.0000 0.0175 0.010802875
## JPM JPM -0.0007657541 1.1344930 0.0000 0.9094 0.011629779
## V V -0.0012712918 0.9563165 0.0000 0.8281 0.010358201
## WMT WMT 0.0035798416 0.4941423 0.0000 0.5466 0.007059840
## XOM XOM -0.0016810638 0.9329538 0.0000 0.8785 0.010191471
## MA MA -0.0014104659 1.0864504 0.0000 0.8381 0.011286918
## PG PG 0.0011782718 0.4169185 0.0004 0.8402 0.006508722
## ORCL ORCL 0.0037151505 1.0051208 0.0000 0.6185 0.010706499
## COST COST 0.0107026790 0.7685684 0.0000 0.1008 0.009018315
## JNJ JNJ -0.0046501256 0.5284799 0.0000 0.3983 0.007304894
## HD HD -0.0008454496 0.9927173 0.0000 0.8973 0.010617980
## MRK MRK 0.0026218085 0.4007275 0.0049 0.7165 0.006393174
## NKE NKE -0.0083129083 1.0534178 0.0000 0.3326 0.011051177
## ABBV ABBV 0.0049295568 0.5753174 0.0007 0.5592 0.007639155
## NFLX NFLX -0.0035636554 1.3400821 0.0000 0.8089 0.013096991
## CVX CVX -0.0055413896 1.0936172 0.0000 0.5602 0.011338065
## KO KO -0.0037630125 0.6105680 0.0000 0.4988 0.007890725
## AMD AMD 0.0138085295 1.6726159 0.0000 0.3639 0.015470159
## QCOM QCOM 0.0051306639 1.2687196 0.0000 0.5814 0.012587704
## ADBE ADBE -0.0033152247 1.2965920 0.0000 0.7225 0.012786618
## PEP PEP -0.0015627830 0.5347994 0.0000 0.7350 0.007349994
## CRM CRM -0.0023180603 1.2686408 0.0000 0.8197 0.012587141
market_expected_return
## [1] 0.01066995
# Calculate the actual return of each stock
actual_returns <- colMeans(returns_data)
# Add Actual return to the given Dataframe
betas_df$ActualReturn <- actual_returns
betas_df
## Symbol Alpha Beta PValue_beta PValue_alpha ExpectedReturn
## MSFT MSFT 0.0101192496 0.8862446 0.0000 0.0593 0.009858125
## AAPL AAPL 0.0116759996 1.2245193 0.0000 0.0947 0.012272263
## INTC INTC -0.0166839405 1.0895639 0.0000 0.1566 0.011309138
## GOOG GOOG 0.0082154804 1.0396709 0.0000 0.2541 0.010953070
## AMZN AMZN 0.0003027082 1.1266987 0.0000 0.9734 0.011574154
## META META 0.0037176481 1.1944103 0.0000 0.7839 0.012057386
## UNH UNH 0.0050677177 0.5726595 0.0000 0.4600 0.007620187
## MCD MCD -0.0033893929 0.7192985 0.0000 0.5155 0.008666694
## ETN ETN 0.0135260809 1.0186252 0.0000 0.0175 0.010802875
## JPM JPM -0.0007657541 1.1344930 0.0000 0.9094 0.011629779
## V V -0.0012712918 0.9563165 0.0000 0.8281 0.010358201
## WMT WMT 0.0035798416 0.4941423 0.0000 0.5466 0.007059840
## XOM XOM -0.0016810638 0.9329538 0.0000 0.8785 0.010191471
## MA MA -0.0014104659 1.0864504 0.0000 0.8381 0.011286918
## PG PG 0.0011782718 0.4169185 0.0004 0.8402 0.006508722
## ORCL ORCL 0.0037151505 1.0051208 0.0000 0.6185 0.010706499
## COST COST 0.0107026790 0.7685684 0.0000 0.1008 0.009018315
## JNJ JNJ -0.0046501256 0.5284799 0.0000 0.3983 0.007304894
## HD HD -0.0008454496 0.9927173 0.0000 0.8973 0.010617980
## MRK MRK 0.0026218085 0.4007275 0.0049 0.7165 0.006393174
## NKE NKE -0.0083129083 1.0534178 0.0000 0.3326 0.011051177
## ABBV ABBV 0.0049295568 0.5753174 0.0007 0.5592 0.007639155
## NFLX NFLX -0.0035636554 1.3400821 0.0000 0.8089 0.013096991
## CVX CVX -0.0055413896 1.0936172 0.0000 0.5602 0.011338065
## KO KO -0.0037630125 0.6105680 0.0000 0.4988 0.007890725
## AMD AMD 0.0138085295 1.6726159 0.0000 0.3639 0.015470159
## QCOM QCOM 0.0051306639 1.2687196 0.0000 0.5814 0.012587704
## ADBE ADBE -0.0033152247 1.2965920 0.0000 0.7225 0.012786618
## PEP PEP -0.0015627830 0.5347994 0.0000 0.7350 0.007349994
## CRM CRM -0.0023180603 1.2686408 0.0000 0.8197 0.012587141
## ActualReturn
## MSFT 0.019977375
## AAPL 0.023948263
## INTC -0.005374803
## GOOG 0.019168550
## AMZN 0.011876863
## META 0.015775035
## UNH 0.012687905
## MCD 0.005277301
## ETN 0.024328956
## JPM 0.010864025
## V 0.009086909
## WMT 0.010639681
## XOM 0.008510407
## MA 0.009876452
## PG 0.007686994
## ORCL 0.014421650
## COST 0.019720994
## JNJ 0.002654768
## HD 0.009772531
## MRK 0.009014982
## NKE 0.002738268
## ABBV 0.012568712
## NFLX 0.009533336
## CVX 0.005796675
## KO 0.004127713
## AMD 0.029278688
## QCOM 0.017718368
## ADBE 0.009471394
## PEP 0.005787211
## CRM 0.010269081
# Runing cross-sectional regression
cross_sectional_model <- lm(ActualReturn ~ Beta, data = betas_df)
summary(cross_sectional_model)
##
## Call:
## lm(formula = ActualReturn ~ Beta, data = betas_df)
##
## Residuals:
## Min 1Q Median 3Q Max
## -0.018238 -0.004282 -0.001747 0.004001 0.012092
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.003249 0.004002 0.812 0.424
## Beta 0.008824 0.004028 2.191 0.037 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.006878 on 28 degrees of freedom
## Multiple R-squared: 0.1463, Adjusted R-squared: 0.1158
## F-statistic: 4.799 on 1 and 28 DF, p-value: 0.03697
plot(betas_df$ActualReturn, betas_df$Beta) +abline(cross_sectional_model)
## integer(0)
par(mfrow=c(1,1))
plot(cross_sectional_model)
## Normality testing: Using Shapiro-Wilk test: The null-hypothesis of this test is that the population is normally distributed.
shapiro.test(cross_sectional_model$residuals) #normal
##
## Shapiro-Wilk normality test
##
## data: cross_sectional_model$residuals
## W = 0.96005, p-value = 0.3108
library(moments)
##
## Attaching package: 'moments'
## The following objects are masked from 'package:PerformanceAnalytics':
##
## kurtosis, skewness
## Normolity Testing: Jaque-bera testing: H_0: normal
jarque.test(cross_sectional_model$residuals)
##
## Jarque-Bera Normality Test
##
## data: cross_sectional_model$residuals
## JB = 0.30952, p-value = 0.8566
## alternative hypothesis: greater
## Testing for Heteroscedasticity
bptest(formula(cross_sectional_model) , data = betas_df , studentize = FALSE) #Because the residuals is normal -> studentize = FALSE
##
## Breusch-Pagan test
##
## data: formula(cross_sectional_model)
## BP = 2.5793, df = 1, p-value = 0.1083
## Autocorrelation and Dynamic Models: Durbin and Watson (1951) (DW) order 1
dwtest(cross_sectional_model)
##
## Durbin-Watson test
##
## data: cross_sectional_model
## DW = 2.4834, p-value = 0.9173
## alternative hypothesis: true autocorrelation is greater than 0
## Autocorrelation and Dynamic Models: Breusch–Godfrey test
bgtest(cross_sectional_model, order=10)
##
## Breusch-Godfrey test for serial correlation of order up to 10
##
## data: cross_sectional_model
## LM test = 7.3742, df = 10, p-value = 0.6897
##Data indicators summary
ind_resid.summary<- c("Average" = mean(cross_sectional_model$residuals) , "Standard Deviation" = sd(cross_sectional_model$residuals), "Skewness" = skewness(cross_sectional_model$residuals), "Kurtosis" = kurtosis(cross_sectional_model$residuals))
ind_resid.summary %>%
kbl(caption = "Table 1: Residuals Summary", escape = TRUE)%>%
kable_classic(full_width = FALSE, html_font = "cambria")
Table 1: Residuals Summary
|
|
x
|
|
Average
|
0.0000000
|
|
Standard Deviation
|
0.0067584
|
|
Skewness
|
-0.2077763
|
|
Kurtosis
|
3.2737391
|