firm_data1 = read.csv('3firmExample_data3.csv')
str(firm_data1)
## 'data.frame': 59 obs. of 4 variables:
## $ date : Factor w/ 59 levels "1995/10/1","1995/11/1",..: 4 5 6 7 8 9 10 1 2 3 ...
## $ Nordstrom: num -0.03615 -0.0568 0.07821 -0.00302 -0.02757 ...
## $ Starbucks: num 0.00521 -0.02105 0.21244 0.2036 0.04797 ...
## $ Microsoft: num 0.1213 0.13923 0.03529 0.06501 0.00138 ...
firm_data1$date
## [1] 1995/3/1 1995/4/1 1995/5/1 1995/6/1 1995/7/1 1995/8/1 1995/9/1
## [8] 1995/10/1 1995/11/1 1995/12/1 1996/1/1 1996/2/1 1996/3/1 1996/4/1
## [15] 1996/5/1 1996/6/1 1996/7/1 1996/8/1 1996/9/1 1996/10/1 1996/11/1
## [22] 1996/12/1 1997/1/1 1997/2/1 1997/3/1 1997/4/1 1997/5/1 1997/6/1
## [29] 1997/7/1 1997/8/1 1997/9/1 1997/10/1 1997/11/1 1997/12/1 1998/1/1
## [36] 1998/2/1 1998/3/1 1998/4/1 1998/5/1 1998/6/1 1998/7/1 1998/8/1
## [43] 1998/9/1 1998/10/1 1998/11/1 1998/12/1 1999/1/1 1999/2/1 1999/3/1
## [50] 1999/4/1 1999/5/1 1999/6/1 1999/7/1 1999/8/1 1999/9/1 1999/10/1
## [57] 1999/11/1 1999/12/1 2000/1/1
## 59 Levels: 1995/10/1 1995/11/1 1995/12/1 1995/3/1 1995/4/1 ... 2000/1/1
library(xts)
## Loading required package: zoo
##
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
##
## as.Date, as.Date.numeric
library(PerformanceAnalytics)
##
## Attaching package: 'PerformanceAnalytics'
## The following object is masked from 'package:graphics':
##
## legend
date1 = as.Date(firm_data1[,1], "%Y/%m/%d")
#convert firm_data1 into time series data: xts
firm_data1.xts = as.xts(firm_data1[,-1], order.by = date1)
firm.data1<-coredata(firm_data1.xts)
summary(firm.data1)
## Nordstrom Starbucks Microsoft
## Min. :-0.212880 Min. :-0.47970 Min. :-0.17634
## 1st Qu.:-0.057395 1st Qu.:-0.01734 1st Qu.:-0.01826
## Median : 0.004950 Median : 0.04064 Median : 0.03848
## Mean : 0.001545 Mean : 0.02846 Mean : 0.04271
## 3rd Qu.: 0.064980 3rd Qu.: 0.09416 3rd Qu.: 0.11174
## Max. : 0.312480 Max. : 0.27967 Max. : 0.28153
skewness(firm.data1)
## Nordstrom Starbucks Microsoft
## Skewness 0.2422885 -0.8882285 0.1712068
rbind(apply(firm.data1, 2, summary),
apply(firm.data1, 2, skewness),
apply(firm.data1, 2, kurtosis))
## Nordstrom Starbucks Microsoft
## Min. -0.212880000 -0.47970000 -0.17634000
## 1st Qu. -0.057395000 -0.01734000 -0.01826500
## Median 0.004950000 0.04064000 0.03848000
## Mean 0.001545085 0.02846068 0.04271153
## 3rd Qu. 0.064980000 0.09415500 0.11173500
## Max. 0.312480000 0.27967000 0.28153000
## 0.242288454 -0.88822851 0.17120676
## 0.351952075 1.85144933 -0.08728437
library(plyr)
library(quantmod)
## Loading required package: TTR
## Version 0.4-0 included new data defaults. See ?getSymbols.
tickers<-c("JWN", "SBUX", "MSFT")
data.env<-new.env()
# here we use l_ply so that we don't double save the data
# getSymbols() does this already so we just want to be memory efficient
# go through every stock and try to use getSymbols()
l_ply(tickers, function(sym) try(getSymbols(sym, env=data.env), silent=T))
## 'getSymbols' currently uses auto.assign=TRUE by default, but will
## use auto.assign=FALSE in 0.5-0. You will still be able to use
## 'loadSymbols' to automatically load data. getOption("getSymbols.env")
## and getOption("getSymbols.auto.assign") will still be checked for
## alternate defaults.
##
## This message is shown once per session and may be disabled by setting
## options("getSymbols.warning4.0"=FALSE). See ?getSymbols for details.
##
## WARNING: There have been significant changes to Yahoo Finance data.
## Please see the Warning section of '?getSymbols.yahoo' for details.
##
## This message is shown once per session and may be disabled by setting
## options("getSymbols.yahoo.warning"=FALSE).
# now we only want the stocks that got stored from getSymbols()
# basically we drop all "bad" tickers
stocks <- tickers[tickers %in% ls(data.env)]
# now we just loop through and merge our good stocks
# if you prefer to use an lapply version here, that is also fine
# since now we are just collecting all the good stock xts() objects
data <- xts()
# i=1
for(i in seq_along(stocks)) {
symbol <- stocks[i]
data <- merge(data, Ad(get(symbol, envir=data.env)))
}
head(data)
## JWN.Adjusted SBUX.Adjusted MSFT.Adjusted
## 2007-01-03 08:00:00 35.71681 14.13297 22.47883
## 2007-01-04 08:00:00 36.32145 14.14900 22.44119
## 2007-01-05 08:00:00 35.64033 14.08886 22.31321
## 2007-01-08 08:00:00 35.82799 14.03674 22.53152
## 2007-01-09 08:00:00 36.21721 13.97660 22.55411
## 2007-01-10 08:00:00 36.39095 13.93249 22.32826
str(data)
## An 'xts' object on 2007-01-03 08:00:00/2019-03-27 08:00:00 containing:
## Data: num [1:3079, 1:3] 35.7 36.3 35.6 35.8 36.2 ...
## - attr(*, "dimnames")=List of 2
## ..$ : NULL
## ..$ : chr [1:3] "JWN.Adjusted" "SBUX.Adjusted" "MSFT.Adjusted"
## Indexed by objects of class: [POSIXct,POSIXt] TZ:
## xts Attributes:
## NULL
# convert POSIXct into date series
data<-xts(coredata(data), order.by = as.Date(index(data), tz=""))
head(data)
## JWN.Adjusted SBUX.Adjusted MSFT.Adjusted
## 2007-01-03 35.71681 14.13297 22.47883
## 2007-01-04 36.32145 14.14900 22.44119
## 2007-01-05 35.64033 14.08886 22.31321
## 2007-01-08 35.82799 14.03674 22.53152
## 2007-01-09 36.21721 13.97660 22.55411
## 2007-01-10 36.39095 13.93249 22.32826
tail(data)
## JWN.Adjusted SBUX.Adjusted MSFT.Adjusted
## 2019-03-20 43.45 71.63 117.52
## 2019-03-21 43.76 72.26 120.22
## 2019-03-22 42.93 71.96 117.05
## 2019-03-25 43.54 72.30 117.66
## 2019-03-26 43.73 72.96 117.91
## 2019-03-27 44.65 72.74 116.77
library(fBasics)
## Loading required package: timeDate
##
## Attaching package: 'timeDate'
## The following objects are masked from 'package:PerformanceAnalytics':
##
## kurtosis, skewness
## Loading required package: timeSeries
##
## Attaching package: 'timeSeries'
## The following object is masked from 'package:zoo':
##
## time<-
##
## Attaching package: 'fBasics'
## The following object is masked from 'package:TTR':
##
## volatility
Sigma = cov(firm_data1[,2:4])
std = sqrt(diag(Sigma))
ones = rep(1,3)
one.vec = matrix(ones, ncol=1)
a = inv(Sigma)%*%one.vec
b = t(one.vec)%*%a
mvp.w =a / as.numeric(b)
mvp.w
## [,1]
## Nordstrom 0.3635998
## Starbucks 0.1936537
## Microsoft 0.4427465
mvp.ret<-sum((mvp.w)*colMeans(firm_data1[,2:4]))
mvp.ret
## [1] 0.02498369
mu<-0.06/12
return <- firm_data1[,2:4]
Ax <- rbind(2*cov(return), colMeans(return), rep(1, ncol(return)))
Ax <- cbind(Ax, rbind(t(tail(Ax, 2)), matrix(0, 2, 2)))
b0 <- c(rep(0, ncol(return)), mu, 1)
out<-solve(Ax, b0)
wgt<-out[1:3]
wgt
## Nordstrom Starbucks Microsoft
## 0.875635380 0.116816458 0.007548163
sum(wgt)
## [1] 1
ret.out<-sum(wgt*colMeans(return))
ret.out.annual<-ret.out*12
ret.out.annual
## [1] 0.06
std.out<-sqrt(t(wgt)%*%cov(return)%*%wgt)
std.out.annual<-std.out*sqrt(12)
std.out.annual
## [,1]
## [1,] 0.335302
return = firm_data1[,2:4]
#specified portfolio return: mu
mu=0.06/12
minvariance <- function(return, mu) {
#return <- log(tail(assets, -1) / head(assets, -1))
Ax <- rbind(2*cov(return), colMeans(return), rep(1, ncol(return)))
Ax <- cbind(Ax, rbind(t(tail(Ax, 2)), matrix(0, 2, 2)))
b0 <- c(rep(0, ncol(return)), mu, 1)
zx<-solve(Ax, b0)
weight<-zx[1:ncol(return)]
ret.out<-sum(weight*colMeans(return))
std.out<-sqrt(t(wgt)%*%cov(return)%*%wgt)
list(weight=weight, rtn=ret.out, sd=std.out)
}
minvariance(return, mu)
## $weight
## Nordstrom Starbucks Microsoft
## 0.875635380 0.116816458 0.007548163
##
## $rtn
## [1] 0.005
##
## $sd
## [,1]
## [1,] 0.09679334
frontier <- function(return){
#return <- log(tail(assets, -1) / head(assets, -1))
n = ncol(return)
Q = cov(return)
Ax <- rbind(2*cov(return), colMeans(return), rep(1, n))
Ax <- cbind(Ax, rbind(t(tail(Ax, 2)), matrix(0, 2, 2)))
r <- colMeans(return)
rbase <- seq(min(r), max(r), length = 100)
s <- sapply(rbase, function(x) {
b0 <- c(rep(0, ncol(return)), x, 1)
y <- head(solve(Ax, b0), n)
sqrt(y%*%Q%*%y)
})
plot(s, rbase, xlab = 'Std', ylab = 'Return')
}
frontier(return)
library(timeSeries)
library(PerformanceAnalytics)