import pandas as pd
import pandas_datareader.data as web
import datetime
import backtrader as bt
import numpy as np
%matplotlib inline
import matplotlib.pyplot as plt
import pyfolio as pf
import quantstats
plt.rcParams["figure.figsize"] = (10, 6) # (w, h)
import sys
from scipy.stats import rankdata
from scipy.stats import stats
from scipy.optimize import minimize

start = '2006-07-21'
end = '2021-05-07'

vti = web.DataReader("VTI", 'yahoo', start, end)['Adj Close'].to_frame("vti")
tlt = web.DataReader("TLT", 'yahoo', start, end)['Adj Close'].to_frame("tlt")
tip = web.DataReader("TIP", 'yahoo', start, end)['Adj Close'].to_frame("tip")
iau = web.DataReader("IAU", 'yahoo', start, end)['Adj Close'].to_frame("iau")
gsg = web.DataReader("GSG", 'yahoo', start, end)['Adj Close'].to_frame("gsg")
shy = web.DataReader("SHY", 'yahoo', start, end)['Adj Close'].to_frame("shy")

price_df = pd.concat([vti, tlt, tip, iau, gsg], axis=1)
return_df = price_df.pct_change().dropna(axis=0)

cash_return_df = shy.pct_change().dropna(axis=0)

return_df.head()

cash_return_df.head()

def RiskParityTargetVol(rets, target, lb, ub) :

    covmat = pd.DataFrame.cov(rets)

    # 목적함수
    def RiskParityObjective(x) :
            
            variance = x.T @ covmat @ x
            sigma = variance ** 0.5
            mrc = 1/sigma * (covmat @ x)
            rc = x * mrc
            a = np.reshape(rc.to_numpy(), (len(rc), 1))
            # 각 자산의 위험 기여도가 같도록
            # 아래쪽의 목적함수가 위험 기여도 차이의 제곱의 합이니 그것이 0에 가깝게 가야 함
            risk_diffs = a - a.T
            # np.ravel: 다차원 배열 -> 1차원 배열
            # 1차원으로 풀어준 후에 제곱 -> 합 -> 목적함수
            sum_risk_diffs_squared = np.sum(np.square(np.ravel(risk_diffs)))
            return (sum_risk_diffs_squared)
    
    # 제약조건: 연간 변동성이 target vol의 95% 이상일 것
    def TargetVolLowerBound(x) :
        
        variance = x.T @ covmat @ x
        sigma = variance ** 0.5 
        sigma_scale = sigma * np.sqrt(252) 
        
        vol_diffs = sigma_scale - (target * 0.95)
        return(vol_diffs)  

    # 제약조건: 연간 변동성이 target vol의 105% 이하일 것
    def TargetVolUpperBound(x) :
        
        variance = x.T @ covmat @ x
        sigma = variance ** 0.5 
        sigma_scale = sigma * np.sqrt(252) 
        
        vol_diffs = (target * 1.05) - sigma_scale
        return(vol_diffs)      
    
    x0 = np.repeat(1/covmat.shape[1], covmat.shape[1]) 
    lbound  = np.repeat(lb, covmat.shape[1])
    ubound  = np.repeat(ub, covmat.shape[1])
    bnds = tuple(zip(lbound, ubound))
    constraints = ({'type': 'ineq', 'fun': TargetVolLowerBound},
                   {'type': 'ineq', 'fun': TargetVolUpperBound})
    options = {'ftol': 1e-20, 'maxiter': 2000}
    
    result = minimize(fun = RiskParityObjective,
                      x0 = x0,
                      method = 'SLSQP',
                      constraints = constraints,
                      options = options,
                      bounds = bnds)
    return(result.x)

# 62 days(3 months) 기준으로 rolling 
weight_df = pd.DataFrame(columns=['vti','tlt','tip','iau','gsg'])

for i in range(len(return_df)-62):
    
    #print(i)
    
    # 63일 분량 데이터 잘라오기 (for문 돌 때마다 index 0~62, 1~63, 2~64, ...)
    temp_return_df = return_df.iloc[i:(i+63), :]
    
    # 63일 중 62일치를 계산에 사용하고 63번째 행은 날짜를 추출
    date_index = str(temp_return_df.iloc[62, :].name)[0:10]
    # Target vol 10% (0.1), 최소/최대 보유 비중은 0과 10 (Long Only, 10배 레버리지까지 허용)
    temp = pd.DataFrame(RiskParityTargetVol(temp_return_df.iloc[0:62, :], 0.1, 0, 10), index=['vti','tlt','tip','iau','gsg'], columns=[date_index]).T
    
    weight_df = pd.concat([weight_df, temp])

weight_df

weight_df_final = weight_df

weight_df_final['shy'] = 1 - weight_df_final.sum(axis = 1)
weight_df_final

weight_df_final.plot(title = 'Weights for each asset class')

sum_of_asset_weight = 1 - weight_df_final['shy']
sum_of_asset_weight.plot(title = 'Sum of weights for each asset group')

rp_tvol_return = (weight_df_final['vti'] * return_df.iloc[62:, :]['vti'] +
                  weight_df_final['tlt'] * return_df.iloc[62:, :]['tlt'] + 
                  weight_df_final['tip'] * return_df.iloc[62:, :]['tip'] +
                  weight_df_final['iau'] * return_df.iloc[62:, :]['iau'] +
                  weight_df_final['gsg'] * return_df.iloc[62:, :]['gsg'] + 
                  weight_df_final['shy'] * cash_return_df.iloc[62:, :]['shy'])
        
rp_tvol_return

rp_tvol_return.index

rp_tvol_return.index = pd.to_datetime(rp_tvol_return.index)

rp_tvol_return.index

quantstats.reports.plots(rp_tvol_return, mode='basic')

quantstats.reports.metrics(rp_tvol_return, mode='full')

quantstats.reports.html(rp_tvol_return, output='Report_RPTVOL.html', title='Risk Parity with 10% Target Vol')