Building optimal portfolios for all combinations of stocks

Question

Hi I recently wrote some code in python that does the following:

1.) Pulls stock closing data from yahoo finance for x number of stock

2.) finds all possible combinations of x stocks in groups of y size (so all combinations of 13 stocks in groups of 10)

3.) applies some calculations to each of these groups

4.) returns the optimal weights of each group. (weights means what percentage of money is to be placed in each stock and must = 100%)

5.) Creates a Dataframe with all the optimal portfolios (groups) weights

My code works however it is slow, each loop on average on my PC takes 2.45 seconds. This is fine for a small number of permutations such as the example above however when the number of selections increase the number of possibilities also increases significantly. For example a list of 30 stocks taken in unique groups of 15 has 155117520 possibilities which would take my code over 12 years..... Just looking for any suggestions or directions to improve the execution speed of my code. I am relatively new to coding but am aware python is slower than other languages at this task however I currently only know some python basics. I use a few for loops in this code which I know are slow and am exploring using .apply() instead, if you could help in any way it would be appreciated.

import pandas as pd
from pandas_datareader import data
import datetime
import numpy as np
import random
import itertools
import requests
import time

time1 = time.time()


start = datetime.datetime(2015, 1, 1)

end = datetime.datetime(2019, 12, 31)

list2 = []

num = []


#################SP Download#######################

sptickers1 = ['MMM', 'ABT', 'ABBV', 'ABMD', 'ACN', 'ATVI', 'ADBE', 'AMD', 'AAP', 'AES', 'AFL', 'A', 'APD', 'AKAM', 'ALK', 'ALB', 'ARE', 'ALXN', 'ALGN', 'ALLE', 'AGN', 'ADS', 'LNT', 'ALL', 'GOOGL', 'GOOG', 'MO', 'AMZN', 'AMCR', 'AEE', 'AAL', 'AEP', 'AXP', 'AIG', 'T', 'AMT', 'AWK', 'AMP', 'ABC', 'AME', 'AMGN', 'APH', 'ADI', 'ANSS', 'ANTM', 'AON', 'AOS', 'APA', 'AIV', 'AAPL', 'AMAT', 'APTV', 'ADM', 'ARNC', 'ANET', 'AJG', 'AIZ', 'ATO', 'ADSK', 'ADP', 'AZO', 'AVB', 'AVY', 'BKR', 'BLL', 'BAC', 'BK', 'BAX', 'BDX', 'BBY', 'BIIB', 'BLK', 'BA', 'BKNG', 'BWA', 'BXP', 'BSX', 'BMY', 'AVGO', 'BR', 'CHRW', 'COG', 'CDNS', 'CPB', 'COF', 'CPRI', 'CAH', 'KMX', 'CCL', 'CAT', 'CBOE', 'CBRE', 'CDW', 'CE', 'CNC', 'CNP', 'CTL', 'CERN', 'CF', 'SCHW', 'CHTR', 'CVX', 'CMG', 'CB', 'CHD', 'CI', 'CINF', 'CTAS', 'CSCO', 'C', 'CFG', 'CTXS', 'CLX', 'CME', 'CMS', 'KO', 'CTSH', 'CL', 'CMCSA', 'CMA', 'CAG', 'CXO', 'COP', 'ED', 'STZ', 'COO', 'CPRT', 'GLW', 'CTVA', 'COST', 'COTY', 'CCI', 'CSX', 'CMI', 'CVS', 'DHI', 'DHR', 'DRI', 'DVA', 'DE', 'DAL', 'XRAY', 'DVN', 'FANG', 'DLR', 'DFS', 'DISCA', 'DISCK', 'DISH', 'DG', 'DLTR', 'D', 'DOV', 'DOW', 'DTE', 'DUK', 'DRE', 'DD', 'DXC', 'ETFC', 'EMN', 'ETN', 'EBAY', 'ECL', 'EIX', 'EW', 'EA', 'EMR', 'ETR', 'EOG', 'EFX', 'EQIX', 'EQR', 'ESS', 'EL', 'EVRG', 'ES', 'RE', 'EXC', 'EXPE', 'EXPD', 'EXR', 'XOM', 'FFIV', 'FB', 'FAST', 'FRT', 'FDX', 'FIS', 'FITB', 'FE', 'FRC', 'FISV', 'FLT', 'FLIR', 'FLS', 'FMC', 'F', 'FTNT', 'FTV', 'FBHS', 'FOXA', 'FOX', 'BEN', 'FCX', 'GPS', 'GRMN', 'IT', 'GD', 'GE', 'GIS', 'GM', 'GPC', 'GILD', 'GL', 'GPN', 'GS', 'GWW', 'HRB', 'HAL', 'HBI', 'HOG', 'HIG', 'HAS', 'HCA', 'PEAK', 'HP', 'HSIC', 'HSY', 'HES', 'HPE', 'HLT', 'HFC', 'HOLX', 'HD', 'HON', 'HRL', 'HST', 'HPQ', 'HUM', 'HBAN', 'HII', 'IEX', 'IDXX', 'INFO', 'ITW', 'ILMN', 'INCY', 'IR', 'INTC', 'ICE', 'IBM', 'IP', 'IPG', 'IFF', 'INTU', 'ISRG', 'IVZ', 'IPGP', 'IQV', 'IRM', 'JKHY', 'J', 'JBHT', 'SJM', 'JNJ', 'JCI', 'JPM', 'JNPR', 'KSU', 'K', 'KEY', 'KEYS', 'KMB', 'KIM', 'KMI', 'KLAC', 'KSS', 'KHC', 'KR', 'LB', 'LHX', 'LH', 'LRCX', 'LW', 'LVS', 'LEG', 'LDOS', 'LEN', 'LLY', 'LNC', 'LIN', 'LYV', 'LKQ', 'LMT', 'L', 'LOW', 'LYB', 'MTB', 'M', 'MRO', 'MPC', 'MKTX', 'MAR', 'MMC', 'MLM', 'MAS', 'MA', 'MKC', 'MXIM', 'MCD', 'MCK', 'MDT', 'MRK', 'MET', 'MTD', 'MGM', 'MCHP', 'MU', 'MSFT', 'MAA', 'MHK', 'TAP', 'MDLZ', 'MNST', 'MCO', 'MS', 'MOS', 'MSI', 'MSCI', 'MYL', 'NDAQ', 'NOV', 'NTAP', 'NFLX', 'NWL', 'NEM', 'NWSA', 'NWS', 'NEE', 'NLSN', 'NKE', 'NI', 'NBL', 'JWN', 'NSC', 'NTRS', 'NOC', 'NLOK', 'NCLH', 'NRG', 'NUE', 'NVDA', 'NVR', 'ORLY', 'OXY', 'ODFL', 'OMC', 'OKE', 'ORCL', 'PCAR', 'PKG', 'PH', 'PAYX', 'PAYC', 'PYPL', 'PNR', 'PBCT', 'PEP', 'PKI', 'PRGO', 'PFE', 'PM', 'PSX', 'PNW', 'PXD', 'PNC', 'PPG', 'PPL', 'PFG', 'PG', 'PGR', 'PLD', 'PRU', 'PEG', 'PSA', 'PHM', 'PVH', 'QRVO', 'PWR', 'QCOM', 'DGX', 'RL', 'RJF', 'RTN', 'O', 'REG', 'REGN', 'RF', 'RSG', 'RMD', 'RHI', 'ROK', 'ROL', 'ROP', 'ROST', 'RCL', 'SPGI', 'CRM', 'SBAC', 'SLB', 'STX', 'SEE', 'SRE', 'NOW', 'SHW', 'SPG', 'SWKS', 'SLG', 'SNA', 'SO', 'LUV', 'SWK', 'SBUX', 'STT', 'STE', 'SYK', 'SIVB', 'SYF', 'SNPS', 'SYY', 'TMUS', 'TROW', 'TTWO', 'TPR', 'TGT', 'TEL', 'FTI', 'TFX', 'TXN', 'TXT', 'TMO', 'TIF', 'TJX', 'TSCO', 'TT', 'TDG', 'TRV', 'TFC', 'TWTR', 'TSN', 'UDR', 'ULTA', 'USB', 'UAA', 'UA', 'UNP', 'UAL', 'UNH', 'UPS', 'URI', 'UTX', 'UHS', 'UNM', 'VFC', 'VLO', 'VAR', 'VTR', 'VRSN', 'VRSK', 'VZ', 'VRTX', 'V', 'VNO', 'VMC', 'WRB', 'WAB', 'WMT', 'WBA', 'DIS', 'WM', 'WAT', 'WEC', 'WFC', 'WELL', 'WDC', 'WU', 'WRK', 'WY', 'WHR', 'WMB', 'WLTW', 'WYNN', 'XEL', 'XRX', 'XLNX', 'XYL', 'YUM', 'ZBRA', 'ZBH', 'ZION', 'ZTS']

dstocks = sptickers1[0:13]

df = data.DataReader(dstocks, 'yahoo', start, end)['Close']

combinations = list(itertools.combinations(dstocks, 10))

combinationslist = []

for i in combinations:
    combinationslist.append(list(i))

for i in combinationslist:
    try:

        start_time = time.time()

        df1 = df[i].copy()

        dfpct = df1.pct_change().apply(lambda x: np.log(x+1))

        sdd = dfpct.std()

        sda = sdd.apply(lambda x: x*np.sqrt(250))

        var = dfpct.var()




        cov_matrix = dfpct.cov()


        dfer = df1.resample('Y').last().pct_change()

        er = dfer.mean()



        p_ret = []
        p_vol = []
        p_weights = []
        num_p = 1000



        for portfolio in range(num_p):
            n = len(i)
            weights = [random.random() for e in range(n)]
            sum_weights = sum(weights)
            weights = [w/sum_weights for w in weights]
            p_weights.append(weights)
            returns = np.dot(weights, er)
            p_ret.append(returns)
            p_var = cov_matrix.mul(weights, axis = 0).mul(weights, axis=1).sum().sum()
            p_sd = np.sqrt(p_var)
            p_sda = p_sd*np.sqrt(250)
            p_vol.append(p_sda)



        data = {'Returns':p_ret, 'Volatility':p_vol}



        for counter, symbol in enumerate(dfpct.columns.tolist()):
            data[symbol+ ' Weight'] = [w[counter] for w in p_weights]



        portfolios = pd.DataFrame(data)



        rf = 0.02

        optimaln = ((portfolios['Returns']-rf)/portfolios['Volatility']).idxmax()

        optimal = portfolios.loc[optimaln]

        optimal1 = pd.DataFrame(optimal).transpose()


        optimal1I = optimal1.index.tolist()

        dictoptimal = portfolios.loc[optimal1I].to_dict(orient='records')

        list2.append(dictoptimal)

        end_time = time.time()

        print("total time taken this loop: ", end_time - start_time)


    except:
        print('Didnt work')
        num.append('didnt work')
        continue

print(len(num))




fin = pd.DataFrame.from_dict(list2)



time2 = time.time()


print('program took ' + str(time2-time1) + ' Seconds')

Answer 1

Checking over 100 million possibilities is going to be slow in any language. That being said, here are some ways to speed up the code a bit.

• There is no need to actually get the list of all combinations when you only need them one at a time. Especially when you have 100 million of them this will be very big in memory. Instead just use it as the generator itertools.combinations returns, you can simply iterate over it.

• There is no need to copy the dataframe, you don't modify it anyway.

• The most important thing is using vectorized functions wherever possible. numpy functions by default work on arrays. So instead of

  dfpct = df1.pct_change().apply(lambda x: np.log(x+1))
  

  use

  dfpct = np.log1p(df1.pct_change())
  

  The numpy.log1p function is the same as numpy.log(1 + x), but more accurate if x is close to zero.

  Similarly, for the random weights:

  weights = np.random.rand(n)
  weights /=  weights.sum()
  

• numpy.sum by default sums along all axis, so doing two sum in a row without specifying an axis is meaningless.

• Don't use a bare except clause. This includes e.g. the user pressing Ctrl+C to abort the program (a very real possibility if your program is going to run 12 years), meaning they have to press that 150 million times. At least use except Exception and then you may as well print the error at least, so you know what is going wrong, by doing except Exception as e. You should also constrain the range of the try..except block as far as possible, e.g. only around the lines you know can cause problems. This is so you don't ignore unexpected errors.

• Using 'Didnt work' as a special value for when something went wrong is maybe not the best idea. Consider using None or np.nan instead.

• Python has an official style-guide, PEP8. It recommends using lower_case for variables and functions and not to use unnecessary whitespace (of which you have a lot).

With these mostly fixed, your code would look like this:

import pandas as pd
from pandas_datareader import data
import datetime
import numpy as np
import random
import itertools
import requests
import time

time1 = time.time()
start = datetime.datetime(2015, 1, 1)
end = datetime.datetime(2019, 12, 31)
list2 = []
num = []

#################SP Download#######################

sptickers1 = ['MMM', 'ABT', 'ABBV', 'ABMD', 'ACN', 'ATVI', 'ADBE', 'AMD', 'AAP', 'AES', 'AFL', 'A', 'APD', ...]
dstocks = sptickers1[:13]
df = data.DataReader(dstocks, 'yahoo', start, end)['Close']
combinations = itertools.combinations(dstocks, 10)

for i in combinations:
    try:
        start_time = time.time()
        df1 = df[list(i)]
        dfpct = np.log1p(df1.pct_change())
        sdd = dfpct.std()
        sda = sdd * np.sqrt(250)
        var = dfpct.var()
        cov_matrix = dfpct.cov()
        dfer = df1.resample('Y').last().pct_change()
        er = dfer.mean()

        p_ret = []
        p_vol = []
        p_weights = []
        num_p = 1000
        for portfolio in range(num_p):
            n = len(i)
            weights = np.random.rand(n)
            weights /=  weights.sum()
            p_weights.append(weights)
            returns = np.dot(weights, er)
            p_ret.append(returns)
            p_var = cov_matrix.mul(weights, axis=0).mul(weights, axis=1).sum()
            p_sda = np.sqrt(p_var)*np.sqrt(250)
            p_vol.append(p_sda)

        data = {'Returns': p_ret, 'Volatility': p_vol}
        for counter, symbol in enumerate(dfpct.columns.tolist()):
            data[symbol+ ' Weight'] = [w[counter] for w in p_weights]

        portfolios = pd.DataFrame(data)
        rf = 0.02
        optimaln = ((portfolios['Returns']-rf)/portfolios['Volatility']).idxmax()
        optimal = portfolios.loc[optimaln].T 
        optimal1I = optimal.index.tolist()
        dictoptimal = portfolios.loc[optimal1I].to_dict(orient='records')
        list2.append(dictoptimal)

        end_time = time.time()
        print("total time taken this loop: ", end_time - start_time)
    except Exception as e:
        print('Didnt work', e)
        num.append('didnt work')
        continue

print(len(num))
fin = pd.DataFrame.from_dict(list2)
time2 = time.time()
print('program took ' + str(time2-time1) + ' Seconds')

This can probably be improved further, but at this point it is hard to follow what exactly your code does. In order to improve this, put independent things into their own function, which allows you to give them a clear name and docstring, explaining what the function does. You should probably have at least a optimal_portfolio and a random_weights function.

Together with this you should also try to come up with more meaningful names than sdd and sda or optimal, optimal1, optimal1I and dictoptimal. Naming things is hard, though.

Once you have done that, you should profile your code to determine which function takes the longest to identify where you need to focus your attention next. The easiest way to do that is to run your script as python -m cProfile -s cumtime script.py.