from math import ceil,floor,isnan
from datetime import datetime
import pandas as pd
import numpy as np
from scipy.optimize import minimize
#from Risk.MaximumDrawdownPercentPerSecurity import MaximumDrawdownPercentPerSecurity
# Based on the following discussion - https://www.quantconnect.com/forum/discussion/3216/classical-asset-allocation-with-mean-variance-optimization/p1
# add data normalization
# no SL and TP as its monthly rebalancing  

# Method
# Each month we estimate the optimal mix of assets weights based on information from the prior 252 trading days and use that mix for the next month. For the covariance matrix, we used the historical covariance matrix of returns for the trailing twelve months. 
# As the mean-variance optimization seeks any optimal set of portfolio weights, There is the potential for the portfolio to become quite concentrated at times. To reduce this possibility, we imposed caps (max weights) on assets to enforce greater diversification as indicated in the paper. E.g. impose a cap of 25% for all risky assets and no cap (i.e. a cap of 100%) for all cash-like assets.
# Optimization
#  For optimization, we tried the different methods like maximizing the Sharpe Ratio, maximizing the return given target volatility stays unchanged (We use 15% as the target volatility). 



class AssetAllocationAlgorithm(QCAlgorithm):
    
    def Initialize(self):
        # frame 1
        # self.SetStartDate(2016, 1, 1)
        # self.SetEndDate(2020, 12, 31)
        # #frame 2
        # self.SetStartDate(2011, 1, 1)
        # self.SetEndDate(2015, 12, 31)
        # # #frame 3
        # self.SetStartDate(2006, 1, 1)
        # self.SetEndDate(2010, 12, 31)
        # #frame 4
        # self.SetStartDate(2001, 1, 1)
        # self.SetEndDate(2005, 12, 31)
        #11 year test
        self.SetStartDate(2010, 1, 1)
        self.SetEndDate(2020, 12, 31)
        self.SetCash(100000)             # Set Strategy Cash
        self.SetSecurityInitializer(lambda x: x.SetDataNormalizationMode(DataNormalizationMode.TotalReturn))
        
        #note without data normalization - it will be factoring splits and dividends 
        # DataNormalizationMode.Adjusted //Factoring in splits and dividends, default setting if no normalization
        # .SplitAdjusted // Just factoring splits, paying dividends as cash
        # .TotalReturn //Adding dividends to asset price
        # .Raw // Price as raw, dividends paid as cash, quantity adjusted on splits


        # tickers = ["SPY", "TLT", "BND", "VGT"] #yl 
        # tickers = ["SPY","QQQ","DIA","IWM","TLT","GLD"] #JW 156

        # tickers = ["SPY", "IVV"] #103.75
        # tickers = ["SPY", "VTI"] #139.96
        # tickers = ["SPY", "VOO"] #125.95
        # tickers = ["SPY", "QQQ"] #131.27
        # tickers = ["SPY", "MLPY"] #-1.4
        # tickers = ["SPY", "VEA"] #50.79
        # tickers = ["SPY", "IEFA"] #47.95
        #tickers = ["SPY", "AGG"] #111.01
        #tickers = ["SPY", "VWO"] #105.34
        # tickers = ["SPY", "IEMG"] #86.49
        # tickers = ["SPY", "VTV"] #61.29
        # tickers = ["SPY", "VUG"] #100.23
        # tickers = ["SPY", "BND"] #108.17
        # tickers = ["SPY", "IJR"] #105.19
        # tickers = ["SPY", "IWM"] #105.26
        # tickers = ["SPY", "IWF"] #105.09
        # tickers = ["SPY", "IJH"] #85.16
        # tickers = ["SPY", "GLD"] #124.93
        #tickers = ["SPY", "VIG"] #119.13
        # tickers = ["SPY", "VCIT"] #99.57
        # tickers = ["SPY", "LQD"] #130.17
        # tickers = ["SPY", "BNDX"] #92.88
        # tickers = ["SPY", "ITOT"] #104.45
        # tickers = ["SPY", "VCSH"] #92.71
        # tickers = ["SPY", "VYM"] #52.98
        # tickers = ["SPY", "BSV"] #110.92 
        # tickers = ["SPY", "USMV"] #132.84
        # tickers = ["SPY", "IAU"] #119.86
        # tickers = ["SPY", "IXUS"] #55.54
        # tickers = ["SPY", "TIP"] #124.19
        # tickers = ["SPY", "MBB"] #82.1
        # tickers = ["SPY", "IGSB"] #93.12
        # tickers = ["SPY", "HYG"] #80.12
        # tickers = ["SPY", "SDY"] #82.85
        # tickers = ["SPY", "SHY"]  #92.29 
        # tickers = ["SPY", "SCHP"] #103.11
        # tickers = ["SPY", "SLV"] #90.70
        # tickers = ["SPY", "SHV"] #64.06 
        # tickers = ["SPY", "IEF"] #113.28
        # tickers = ["SPY", "TLT"] #131.87

        # tickers = ["IEF", "TLT", "SPY", "EFA", "EEM", "JPXN", "VGT"] #original 75.65%
#by return 130 above 
         
        # tickers = ["SPY", "VTI", "QQQ", "USMV", "LQD", "TLT", "GLD", "IAU"] #scenario1 136.57
        # tickers = ["SPY", "VTI", "USMV", "TLT", "QQQ", "LQD"] scenario3  112.67?
        # tickers = ["SPY", "QQQ", "TLT", "LQD", "GLD"] scenario4 145.33?
        # tickers = ["SPY", "QQQ", "LQD", "GLD", "IAU"] scenario5 129.11?
        # tickers = ["SPY", "VTI", "TLT", "GLD"] #scenario 6 134.61
        # tickers = ["SPY", "USMV", "TLT", "GLD"]  scenario 7 92.98?
        # tickers = ["SPY", "QQQ", "TLT", "IAU"] scenario 10 137?

    
        # tickers = ["SPY", "VTI", "QQQ",  "LQD", "TLT", "GLD", "IAU"] #scenario2 206.90  v7% 122.53 v10% 154.36, v15% 188.23
        # tickers = ["SPY", "QQQ", "TLT", "GLD"] #scenario 8 179.02 v7% 139.84 v10% 165.46, v15% 176.20
        tickers = ["SPY", "QQQ", "LQD", "GLD"] #scenario 9 177.81 v7% 160.86 v10% 159.09 v15% 144.35
        

        
        self.symbols = [] 
        for i in tickers:
            self.symbols.append(self.AddEquity(i, Resolution.Daily).Symbol)
        for syl in self.symbols:
            syl.window = RollingWindow[TradeBar](252)
            # self.SetWarmup(252)
        
        self.Schedule.On(self.DateRules.MonthStart("SPY"), self.TimeRules.AfterMarketOpen("SPY"), Action(self.Rebalancing))
        
        
    ## Risk model - additional parameter
        # stopRisk = self.GetParameter("stopRisk")
        # if stopRisk is None:
        #     stopRisk = 0.16 
        # self.SetRiskManagement(TrailingStopRiskManagementModel(float(stopRisk)))


    def OnData(self, data):
        if data.ContainsKey("SPY"):
            for syl in self.symbols:
                syl.window.Add(data[syl])
        
    def Rebalancing(self):
        
        data = {}
        for syl in self.symbols:
            data[syl] = [float(i.Close) for i in syl.window]
        df_price = pd.DataFrame(data,columns=data.keys()) 
        daily_return = (df_price / df_price.shift(1)).dropna()
            
        a = PortfolioOptimization(daily_return, 0, len(data))
        opt_weight = a.opt_portfolio()  
            
        if isnan(sum(opt_weight)): return
        self.Log(str(opt_weight))
        
        for i in range(len(data)):
            self.SetHoldings(df_price.columns[i], opt_weight[i])
            
            # equally weighted
            # self.SetHoldings(self.symbols[i], 1.0/len(data))


        
class PortfolioOptimization(object):

    import numpy as np
    import pandas as pd

    def __init__(self, df_return, risk_free_rate, num_assets):
        
        self.daily_return = df_return
        self.risk_free_rate = risk_free_rate
        self.n = num_assets # numbers of risk assets in portfolio
        self.target_vol = 0.05   #orig is 0.05 and then we test 0.07, 0.10, 0.15

    def annual_port_return(self, weights):
        # calculate the annual return of portfolio
        return np.sum(self.daily_return.mean() * weights) * 252

    def annual_port_vol(self, weights):
        # calculate the annual volatility of portfolio
        return np.sqrt(np.dot(weights.T, np.dot(self.daily_return.cov() * 252, weights)))

    def min_func(self, weights):
        
        # method 1: maximize sharp ratio
        # return - self.annual_port_return(weights) / self.annual_port_vol(weights)
        
        # # method 2: maximize the return with target volatility 
        return - self.annual_port_return(weights) / self.target_vol

    def opt_portfolio(self):
        # maximize the sharpe ratio to find the optimal weights
        cons = ({'type': 'eq', 'fun': lambda x: np.sum(x) - 1})
        bnds = tuple((0, 1) for x in range(2)) + tuple((0, 0.25) for x in range(self.n - 2))
        opt = minimize(self.min_func,                               # object function
                      np.array(self.n * [1. / self.n]),            # initial value
                      method='SLSQP',                              # optimization method
                      bounds=bnds,                                 # bounds for variables 
                      constraints=cons)                            # constraint conditions
                      
        opt_weights = opt['x']
 
        return opt_weights