SciPy產品組合優化與行業界限分組

Question

嘗試優化投資組合權重分配，通過限制風險最大化我的回報函數。我沒有任何問題可以通過簡單的約束條件找到最優化的權重給我的收益函數，即所有權重之和等於1，並且使我的總風險低於目標風險的其他約束。 我的問題是，如何爲每個組添加行業權重界限？ 我的代碼如下：

# -*- coding: utf-8 -*- 
import pandas as pd 
import numpy as np 
import scipy.optimize as sco 

dates = pd.date_range('1/1/2000', periods=8) 
industry = ['industry', 'industry', 'utility', 'utility', 'consumer'] 
symbols = ['A', 'B', 'C', 'D', 'E'] 
zipped = list(zip(industry, symbols)) 
index = pd.MultiIndex.from_tuples(zipped) 

noa = len(symbols) 

data = np.array([[10, 9, 10, 11, 12, 13, 14, 13], 
       [11, 11, 10, 11, 11, 12, 11, 10], 
       [10, 11, 10, 11, 12, 13, 14, 13], 
       [11, 11, 10, 11, 11, 12, 11, 11], 
       [10, 11, 10, 11, 12, 13, 14, 13]]) 

market_to_market_price = pd.DataFrame(data.T, index=dates, columns=index) 

rets = market_to_market_price/market_to_market_price.shift(1) - 1.0 
rets = rets.dropna(axis=0, how='all') 

expo_factor = np.ones((5,5)) 
factor_covariance = market_to_market_price.cov() 
delta = np.diagflat([0.088024, 0.082614, 0.084237, 0.074648, 
           0.084237]) 
cov_matrix = np.dot(np.dot(expo_factor, factor_covariance), 
          expo_factor.T) + delta 

def calculate_total_risk(weights, cov_matrix): 
    port_var = np.dot(np.dot(weights.T, cov_matrix), weights) 
    return port_var 

def max_func_return(weights): 
    return -np.sum(rets.mean() * weights) 

# optimized return with given risk 
tolerance_risk = 27 
noa = market_to_market_price.shape[1] 
cons = ({'type': 'eq', 'fun': lambda x: np.sum(x) - 1}, 
     {'type': 'eq', 'fun': lambda x: calculate_total_risk(x, cov_matrix) - tolerance_risk}) 
bnds = tuple((0, 1) for x in range(noa)) 
init_guess = noa * [1./noa,] 
opts_mean = sco.minimize(max_func_return, init_guess, method='SLSQP', 
         bounds=bnds, constraints=cons) 


In [88]: rets 
Out[88]: 
      industry    utility   consumer 
        A   B   C   D   E 
2000-01-02 -0.100000 0.000000 0.100000 0.000000 0.100000 
2000-01-03 0.111111 -0.090909 -0.090909 -0.090909 -0.090909 
2000-01-04 0.100000 0.100000 0.100000 0.100000 0.100000 
2000-01-05 0.090909 0.000000 0.090909 0.000000 0.090909 
2000-01-06 0.083333 0.090909 0.083333 0.090909 0.083333 
2000-01-07 0.076923 -0.083333 0.076923 -0.083333 0.076923 
2000-01-08 -0.071429 -0.090909 -0.071429 0.000000 -0.071429 

In[89]: opts_mean['x'].round(3) 
Out[89]: array([ 0.233, 0.117, 0.243, 0.165, 0.243]) 

我怎麼可以添加這樣的基團與5個資產下降，使得總和爲下面的約束？

model = pd.DataFrame(np.array([.08,.12,.05]), index= set(industry), columns = ['strategic']) 
model['tactical'] = [(.05,.41), (.2,.66), (0,.16)] 
In [85]: model 
Out[85]: 
      strategic  tactical 
industry  0.08 (0.05, 0.41) 
consumer  0.12 (0.2, 0.66) 
utility  0.05  (0, 0.16) 

我已閱讀本類似的帖子SciPy optimization with grouped bounds但仍不能得到任何線索，任何機構可以幫助？ 謝謝。

Answer 1

首先，考慮使用專門爲凸優化設計的模塊cvxopt。我不太熟悉，但有效邊界的例子是here。

現在回到你的問題，這是一個解決方法，專門用於你發佈的問題，並使用minimize。 （可以概括爲在輸入類型和用戶友好性方面創建更多的靈活性，並且基於類的實現也將在此處有用。）

關於您的問題，「我如何添加組邊界？」，答案是，你實際上需要通過constraints要做到這一點，而不是bounds參數，因爲

可選的，下限和上限每個元素在X也 使用邊界參數指定。 [着重點]

此規範與您要做的事不符。相反，下面的例子是爲每個組的上限和下限分別添加一個不等式約束。函數mapto_constraints返回添加到當前約束的字典列表。

首先，這裏的一些示例數據：

import pandas as pd 
import numpy as np 
import numpy.random as npr 
npr.seed(123) 
from scipy.optimize import minimize 

# Create a DataFrame of hypothetical returns for 5 stocks across 3 industries, 
# at daily frequency over a year. Note that these will be in decimal 
# rather than numeral form. (i.e. 0.01 denotes a 1% return) 

dates = pd.bdate_range(start='1/1/2000', end='12/31/2000') 
industry = ['industry'] * 2 + ['utility'] * 2 + ['consumer'] 
symbols = list('ABCDE') 
zipped = list(zip(industry, symbols)) 
cols = pd.MultiIndex.from_tuples(zipped) 
returns = pd.DataFrame(npr.randn(len(dates), len(cols)), index=dates, columns=cols) 
returns /= 100 + 3e-3 #drift term 

returns.head() 
Out[191]: 
      industry   utility   consumer 
        A  B  C  D  E 
2000-01-03 -0.01484 0.00986 -0.00476 0.00235 -0.00630 
2000-01-04 0.00518 0.00958 -0.01210 -0.00814 -0.01664 
2000-01-05 0.00233 -0.01665 -0.00366 0.00520 0.02058 
2000-01-06 0.00368 0.01253 0.00259 0.00309 -0.00211 
2000-01-07 -0.00383 0.01174 0.00375 0.00336 -0.00608 

你可以看到年度人物「道理」：

(1 + returns.mean()) ** 252 - 1 
Out[199]: 
industry A -0.05531 
      B 0.32455 
utility C 0.10979 
      D 0.14339 
consumer E -0.12644 

現在對於將在優化中使用的一些功能。這些是從伊夫Hilpisch的Python for Finance例子後藍本，第11章

def logrels(rets): 
    """Log of return relatives, ln(1+r), for a given DataFrame rets.""" 
    return np.log(rets + 1) 

def statistics(weights, rets): 
    """Compute expected portfolio statistics from individual asset returns. 

    Parameters 
    ========== 
    rets : DataFrame 
     Individual asset returns. Use numeral rather than decimal form 
    weights : array-like 
     Individual asset weights, nx1 vector. 

    Returns 
    ======= 
    list of (pret, pvol, pstd); these are *per-period* figures (not annualized) 
     pret : expected portfolio return 
     pvol : expected portfolio variance 
     pstd : expected portfolio standard deviation 

    Note 
    ==== 
    Note that Modern Portfolio Theory (MPT), being a single-period model, 
    works with (optimizes using) continuously compounded returns and 
    volatility, using log return relatives. The difference between these and 
    more commonly used geometric means will be negligible for small returns. 
    """ 

    if isinstance(weights, (tuple, list)): 
     weights = np.array(weights) 
    pret = np.sum(logrels(rets).mean() * weights) 
    pvol = np.dot(weights.T, np.dot(logrels(rets).cov(), weights)) 
    pstd = np.sqrt(pvol) 
    return [pret, pvol, pstd] 

# The below are a few convenience functions around statistics() above, needed 
# because scipy minimize must optimize a function that returns a scalar 

def port_ret(weights, rets): 
    return -1 * statistics(weights=weights, rets=rets)[0] 

def port_variance(weights, rets): 
    return statistics(weights=weights, rets=rets)[1] 

這裏是一個等權重資產組合預期年化標準差。我只是在此處將其用作優化中的錨（risk_tol參數）。

statistics([0.2] * 5, returns)[2] * np.sqrt(252) # ew anlzd stdev 
Out[192]: 0.06642120658640735 

下一個函數需要一個數據幀，看起來像你的model數據幀，並建立約束每個組。請注意，這非常不靈活，因爲您需要遵循特定的返回格式和您現在使用的數據幀。

def mapto_constraints(rets, model): 
    tactical = model['tactical'].to_dict() # values are tuple bounds 
    industries = rets.columns.get_level_values(0) 
    group_cons = list() 
    for key in tactical: 
     if isinstance(industries.get_loc('consumer'), int): 
      pos = [industries.get_loc(key)] 
     else: 
      pos = np.where(industries.get_loc(key))[0].tolist() 
     lb = tactical[key][0] 
     ub = tactical[key][1] # upper and lower bounds 
     lbdict = {'type': 'ineq', 
        'fun': lambda x: np.sum(x[pos[0]:(pos[-1] + 1)]) - lb} 
     ubdict = {'type': 'ineq', 
        'fun': lambda x: ub - np.sum(x[pos[0]:(pos[-1] + 1)])} 
     group_cons.append(lbdict); group_cons.append(ubdict) 
    return group_cons 

的說明，關於如何約束被上述構建的：

等式約束意味着，約束函數結果是成爲零 而不等式意味着它是非負的。

最後，優化自身：

def opt(rets, risk_tol, model, round=3):  
    noa = len(rets.columns) 
    guess = noa * [1./noa,] # equal-weight; needed for initial guess 
    bnds = tuple((0, 1) for x in range(noa)) 
    cons = [{'type': 'eq', 'fun': lambda x: np.sum(x) - 1.}, 
      {'type': 'ineq', 'fun': lambda x: risk_tol - port_variance(x, rets=rets)} 
      ] + mapto_constraints(rets=rets, model=model) 
    opt = minimize(port_ret, guess, args=(returns,), method='SLSQP', bounds=bnds, 
        constraints=cons, tol=1e-10) 
    return opt.x.round(round) 

model = pd.DataFrame(np.array([.08,.12,.05]), 
        index= set(industry), columns = ['strategic']) 
model['tactical'] = [(.05,.41), (.2,.66), (0,.16)] 

# Set variance threshold equal to the equal-weighted variance 
# Note that I set variance as an inequality rather than equality (i.e. 
# resulting variance should be less than threshold). 

opt(returns, risk_tol=port_variance([0.2] * 5, returns), model=model) 
Out[195]: array([ 0.188, 0.225, 0.229, 0.197, 0.16 ])