[Uwgarp-commits] r17 - in pkg/GARPFRM: data sandbox

[noreply at r-forge.r-project.org]

Author: rossbennett34
Date: 2013-11-26 08:02:51 +0100 (Tue, 26 Nov 2013)
New Revision: 17

Added:
   pkg/GARPFRM/data/fama_french_factors.rda
Modified:
   pkg/GARPFRM/sandbox/regression_example.R
Log:
Adding dataset for fama french 3 factors. Adding material to regression example

Added: pkg/GARPFRM/data/fama_french_factors.rda
===================================================================
(Binary files differ)


Property changes on: pkg/GARPFRM/data/fama_french_factors.rda
___________________________________________________________________
Added: svn:mime-type
   + application/octet-stream

Modified: pkg/GARPFRM/sandbox/regression_example.R
===================================================================
--- pkg/GARPFRM/sandbox/regression_example.R	2013-11-26 03:02:08 UTC (rev 16)
+++ pkg/GARPFRM/sandbox/regression_example.R	2013-11-26 07:02:51 UTC (rev 17)
@@ -5,7 +5,110 @@
 # The returns data includes weekly returns
 data(returns)
 
-# Get the weekly returns of AAPL and SPY from the returns object
+# Extract the SPY returns from the returns object
+SPY.ret <- returns[, "SPY"]
+
+# Exploratory Data Analysis, Basic Probability, Basic Statistics
+
+# Plot the SPY weekly returns
+plot(SPY.ret, main="SPY Weekly Returns")
+
+# Plot the density of SPY Weekly Returns
+plot(density(SPY.ret), main="Density of SPY Weekly Returns")
+rug(SPY.ret)
+# sample estimates
+curve(dnorm(x, mean=mean(SPY.ret), sd=sd(SPY.ret)), 
+      add=TRUE, col="red", lty=2, lwd=2)
+# robust estimates
+curve(dnorm(x, mean=median(SPY.ret), sd=mad(SPY.ret)), 
+      add=TRUE, col="blue", lty=2, lwd=2)
+legend("topleft", legend=c("estimated density", "normal density", "robust normal density"), 
+       col=c("black", "red", "blue"), lty=c(1, 2, 2), bty="n", cex=0.8)
+
+# Quantile-Quantile Plot of SPY Weekly Returns
+qqnorm(SPY.ret)
+qqline(SPY.ret)
+
+# Shapiro-Wilk Normality test
+# Null hypothesis that a sample came from a normally distributed population
+shapiro.test(coredata(SPY.ret))
+
+# Sample mean of SPY return
+mean(SPY.ret)
+
+# Sample Variance of SPY returns
+var(SPY.ret)
+
+# Sample standard deviation of SPY returns
+sd(SPY.ret)
+
+# Standard error of SPY returns
+sd(SPY.ret) / sqrt(nrow(SPY.ret))
+
+# Sample skewness of SPY returns
+# See ?skewness for additional methods for calculating skewness
+skewness(SPY.ret, method="sample")
+
+# Sample kurtosis of SPY returns
+# See ?kurtosis for additional methods for calculating kurtosis
+kurtosis(SPY.ret, method="sample")
+
+# Summary statistics of SPY returns
+summary(SPY.ret)
+
+# Sample quantiles of SPY returns
+quantile(SPY.ret, probs=c(0, 0.25, 0.5, 0.75, 1))
+
+# Sample correlation of SPY returns
+cor(returns)
+
+# Sample covariance of SPY returns
+cov(returns)
+
+# Distributions
+# R has functions to compute the density, distribution function, quantile, and 
+# random number generation for several distributions.
+
+# Normal Distribution
+# dnorm, pnorm, qnorm, rnorm
+
+# Chi-Squared Distribution
+# dchisq, pchisq, qchisq, rchisq
+
+# Student t Distribution
+# dt, pt, qt, rt
+
+# F Distribution
+# df, pf, qf, rf
+
+# These functions will be demonstrated for the normal distribution.
+
+# Use dnorm to plot the pdf of a standard normal distribution
+curve(dnorm(x), from=-4, to=4, main="Standard Normal pdf")
+
+# Calculate the probability that Y <= 2 when Y is distributed N(1, 4) with mean of 1 and variance of 4
+pnorm(q=2, mean=1, sd=2)
+pnorm(q=0.5)
+
+# Quantile function of a standard normal at probability 0.975
+qnorm(p=0.975)
+
+# Generate 10 random numbers from a normal distribution with mean 0.0015 and standard deviation 0.025
+rnorm(n=10, mean=0.0015, sd=0.025)
+
+# Hypothesis testing
+# The null hypothesis is that the true mean return of SPY is equal to 0
+t.test(x=SPY.ret, alternative="two.sided", mu=0)
+
+# Replicate the results of t.test using the method outlined in the book
+t_stat <- (mean(SPY.ret) - 0) / (sd(SPY.ret) / sqrt(nrow(SPY.ret)))
+p_value <- 2 * pt(q=-abs(t_stat), df=462)
+df <- nrow(SPY.ret) - 1
+ci <- mean(SPY.ret) + c(-1, 1) * 1.96 * sd(SPY.ret) / sqrt(nrow(SPY.ret))
+
+##### Regression #####
+# Signle Regressor
+# Extract the weekly returns of AAPL and SPY from the returns object
 AAPL.ret <- returns[, "AAPL"]
 SPY.ret <- returns[, "SPY"]
 
@@ -15,7 +118,6 @@
 
 # Fit the linear regression model
 model.fit <- lm(AAPL.ret ~ SPY.ret)
-class(model.fit)
 
 # The print method displays the call and the coefficients of the model
 print(model.fit)
@@ -26,10 +128,8 @@
 resid(model.fit)
 rstandard(model.fit)
 
-# The summary method displays much more information
+# The summary method displays additional information for the linear model
 model.summary <- summary(model.fit)
-class(model.summary)
-names(model.summary)
 
 print(model.summary)
 
@@ -52,65 +152,48 @@
 lines(x=coredata(SPY.ret), y=model.pi[, "upr"], col="red", lty=2)
 lines(x=coredata(SPY.ret), y=model.pi[, "lwr"], col="red", lty=2)
 
+# Plot the residuals of the linear model
 plot(resid(model.fit), type="h")
 
+# Multiple Regressors
+# Fama French Factor Model
+data(returns)
+data(fama_french_factors)
 
-##### EDA #####
+# Fit the factors and AAPL returns
+# AAPL returns can be explained by the Fama French 3-Factor model 
+# AAPL excess returns
 
-# mean returns
-mean(AAPL.ret)
-mean(SPY.ret)
+ff_factors <- fama_french_factors[, 1:3]
 
-# standard deviation
-sd(AAPL.ret)
-sd(SPY.ret)
+# Align the dates of the Fama-French Factors and the returns
+returns <- returns['/2013-10-25']
+returns <- returns[, -1]
+AAPL.ret <- returns[, "AAPL"]
 
-# skewness
-skewness(x=AAPL.ret, method="sample")
-skewness(x=SPY.ret, method="sample")
+AAPL.e <- AAPL.ret - fama_french_factors[, "RF"] / 100
 
-# kurtosis
-kurtosis(x=AAPL.ret, method="sample_excess")
-kurtosis(x=SPY.ret, method="sample_excess")
+ff.fit <- lm(AAPL.e ~ ff_factors)
+print(ff.fit)
+summary(ff.fit)
 
-hist(AAPL.ret)
-hist(SPY.ret)
+# Fit the Fama-French 3 Factor Model to all the assets in the returns object
+# Excess returns
+ret.e <- returns - (fama_french_factors[, "RF"] / 100) %*% rep(1, ncol(returns))
+ff.fit <- lm(ret.e ~ ff_factors)
+print(ff.fit)
+summary(ff.fit)
 
-plot(density(AAPL.ret))
-plot(density(SPY.ret))
 
-mu.est <- mean(AAPL.ret)
-sd.est <- sd(AAPL.ret)
-plot(density(AAPL.ret))
-rug(AAPL.ret)
-# sample estimates
-curve(dnorm(x, mean=mean(AAPL.ret), sd=sd(AAPL.ret)), 
-      add=TRUE, col="red", lty=2, lwd=2)
-# robust estimates
-curve(dnorm(x, mean=median(AAPL.ret), sd=mad(AAPL.ret)), 
-      add=TRUE, col="blue", lty=2, lwd=2)
+beta0 <- coef(ff.fit)[1,]
+beta1 <- coef(ff.fit)[2,]
+beta2 <- coef(ff.fit)[3,]
+beta3 <- coef(ff.fit)[4,]
+rsq <- sapply(X=summary(ff.fit), FUN=function(x) x$r.squared)
+names(rsq) <- colnames(ret.e)
 
-plot(density(SPY.ret))
-rug(SPY.ret)
-# sample estimates
-curve(dnorm(x, mean=mean(SPY.ret), sd=sd(SPY.ret)), 
-      add=TRUE, col="red", lty=2, lwd=2)
-# robust estimates
-curve(dnorm(x, mean=median(SPY.ret), sd=mad(SPY.ret)), 
-      add=TRUE, col="blue", lty=2, lwd=2)
-
-qqnorm(AAPL.ret)
-qqline(AAPL.ret)
-
-boxplot(coredata(returns))
-
-# Shapiro-Wilk Normality test
-# Null hypothesis that a sample came from a normally distributed population
-shapiro.test(coredata(AAPL.ret))
-shapiro.test(coredata(SPY.ret))
-
-#library(MASS)
-#fitdistr(x=AAPL.ret, densfun="normal")
-#fitdistr(x=AAPL.ret, densfun="t", start=list(m=0.001, s=0.05, df=2), lower=c(-0.01, 0, 0.1))
-
-
+par(mfrow=c(2,2))
+barplot(beta1, main="Beta for Market-RF", col=c(2:6))
+barplot(beta2, main="Beta for SMB", col=c(2:6))
+barplot(beta3, main="Beta for HML", col=c(2:6))
+barplot(rsq, main="R Squared Values", col=c(2:6))