a1 <- read.table("AAPL.csv", sep=",", header=TRUE)
a2 <- read.table("C.csv", sep=",", header=TRUE)
a3 <- read.table("CAT.csv", sep=",", header=TRUE)
a4 <- read.table("IBM.csv", sep=",", header=TRUE)
a5 <- read.table("XOM.csv", sep=",", header=TRUE)

r1 <- (a1$Adj.Close[-1]-a1$Adj.Close[-nrow(a1)])/a1$Adj.Close[-nrow(a1)]

r2 <- (a2$Adj.Close[-1]-a2$Adj.Close[-nrow(a2)])/a2$Adj.Close[-nrow(a2)]

r3 <- (a3$Adj.Close[-1]-a3$Adj.Close[-nrow(a3)])/a3$Adj.Close[-nrow(a3)]

r4 <- (a4$Adj.Close[-1]-a4$Adj.Close[-nrow(a4)])/a4$Adj.Close[-nrow(a4)]

r5 <- (a5$Adj.Close[-1]-a5$Adj.Close[-nrow(a5)])/a5$Adj.Close[-nrow(a5)]

ret <- as.data.frame(cbind(r1,r2,r3,r4,r5))



#Use only two stocks:
rr <- ret[,1:2]

#Compute variance-covariance matgrix:
sigma <- cov(rr)

#Compute mean vector:
mm <- colMeans(rr)



#Construct many portfolios:
j <- 0
return_p <- rep(10000)
sd_p <- rep(0,10000)
vect_0 <- rep(0, 10000)
fractions <- matrix(vect_0, 5000,2)


for (a in seq(-1, 1, 0.01)) {
	for (b in seq(-1, 1, 0.01)) {
		if(a+b==1) {
			j=j+1
        fractions[j,] <- c(a,b)
         sd_p[j] <-  (t(fractions[j,]) %*% sigma %*% fractions[j,])^.5
        return_p[j] <-  fractions[j,] %*% mm

			}
			}
			}
						
R_p <- return_p[1:j]
sigma_p <- sd_p[1:j]

plot(sigma_p, R_p, type="l")