Question 1

library(quadprog)

returns_data <- matrix(c(
  0.01, 0.02, -0.01, 0.005,  # Daily returns for ETF "0050" (complete the sequence)
  0.015, 0.018, -0.012, 0.008,  # Daily returns for ETF "0056" (complete the sequence)
  0.012, 0.025, -0.015, 0.007,  # Daily returns for ETF "006205" (complete the sequence)
  0.009, 0.021, -0.008, 0.006  # Daily returns for ETF "00646" (complete the sequence)
), nrow = 4, byrow = TRUE)

colnames(returns_data) <- c("0050", "0056", "006205", "00646")

# Calculate covariance matrix and mean returns
cov_matrix <- cov(returns_data)
mean_returns <- colMeans(returns_data)

# Add small positive value to the diagonal of covariance matrix for regularization
epsilon <- 1e-10
cov_matrix_reg <- cov_matrix + epsilon * diag(ncol(returns_data))

# Define objective function for optimization (minimize portfolio risk)
objective_function <- function(w, cov_matrix) {
  portfolio_risk <- sqrt(t(w) %*% cov_matrix %*% w)
  return(portfolio_risk)
}

# Constraints for optimization (sum of weights equals 1)
constraints <- matrix(1, nrow = 1, ncol = ncol(returns_data))
bounds <- list(lower = rep(0, ncol(returns_data)), upper = rep(1, ncol(returns_data)))

# Perform optimization to find optimal weights
optimal_weights <- solve.QP(Dmat = 2 * cov_matrix_reg, dvec = rep(0, ncol(cov_matrix)),
                            Amat = t(constraints), bvec = 1, meq = 1)$solution

# Calculate return and standard deviation of the GMVP
GMVP_return <- sum(optimal_weights * mean_returns)
GMVP_std_dev <- sqrt(t(optimal_weights) %*% cov_matrix %*% optimal_weights)

print(optimal_weights)
## [1]  0.7336424  0.4253957  0.4679659 -0.6270040
print(GMVP_return)
## [1] 0.008030055
print(GMVP_std_dev)
##              [,1]
## [1,] 2.860615e-07

Question 2

# Sample monthly returns data for the specified period (replace with actual data)
returns_data_monthly <- matrix(c(
  0.05, 0.06, -0.03, 0.04,  # Monthly returns for ETF "0050"
  0.08, 0.07, -0.02, 0.05,  # Monthly returns for ETF "0056"
  0.06, 0.09, -0.04, 0.03,  # Monthly returns for ETF "006205"
  0.07, 0.08, -0.03, 0.06  # Monthly returns for ETF "00646"
), nrow = 4, byrow = TRUE)

colnames(returns_data_monthly) <- c("0050", "0056", "006205", "00646")

# Calculate covariance matrix and mean returns for monthly data
cov_matrix_monthly <- cov(returns_data_monthly)
mean_returns_monthly <- colMeans(returns_data_monthly)

# Add small positive value to the diagonal of covariance matrix for regularization
epsilon <- 1e-10
cov_matrix_reg_monthly <- cov_matrix_monthly + epsilon * diag(ncol(returns_data_monthly))

# Perform optimization to find optimal weights for monthly data
optimal_weights_monthly <- solve.QP(Dmat = 2 * cov_matrix_reg_monthly, dvec = rep(0, ncol(cov_matrix_monthly)),
                                    Amat = t(constraints), bvec = 1, meq = 1)$solution

# Calculate return and standard deviation of the GMVP for monthly data
GMVP_return_monthly <- sum(optimal_weights_monthly * mean_returns_monthly)
GMVP_std_dev_monthly <- sqrt(t(optimal_weights_monthly) %*% cov_matrix_monthly %*% optimal_weights_monthly)

print(optimal_weights_monthly)
## [1] -4.999983e-01  4.999996e-01  9.999981e-01  5.624986e-07
print(GMVP_return_monthly)
## [1] -0.02499984
print(GMVP_std_dev_monthly)
##              [,1]
## [1,] 1.704729e-08

Question 3

# Calculate excess returns (monthly)
rf <- 0  # Risk-free rate
excess_returns <- returns_data_monthly - rf

# Calculate mean excess returns and covariance matrix
mean_excess_returns <- colMeans(excess_returns)
cov_matrix_monthly <- cov(excess_returns)

# Calculate tangency portfolio weights with regularization
tangency_weights <- solve(cov_matrix_monthly + epsilon * diag(ncol(cov_matrix_monthly)), mean_excess_returns) / 
  sum(solve(cov_matrix_monthly + epsilon * diag(ncol(cov_matrix_monthly)), mean_excess_returns))

# Calculate tangency portfolio return and standard deviation
tangency_return <- sum(tangency_weights * mean_excess_returns)
# Calculate tangency portfolio standard deviation
tangency_std_dev <- sqrt(tangency_weights %*% cov_matrix_monthly %*% tangency_weights)

print(tangency_weights)
##          0050          0056        006205         00646 
## -5.000028e-01  4.999994e-01  1.000006e+00 -2.362507e-06
print(tangency_return)
## [1] -0.02500051
print(tangency_std_dev)
##              [,1]
## [1,] 4.854012e-08