Quantopian's community platform is shutting down. Please read this post for more information and download your code.
Back to Community
Help with CVXPY implementing optimization
Shiv Chawla
edited May 24, 2019
Seeking Help Tools and Tips Collaboration Request Portfolio Optimization

Hi everyone

I am trying to create open source version of optimization module provided by the Q team. I have modeled the problem but it is giving me some matrix issue (related to rank). I am not great at python but I am sure that optimization setup is correct (or may be not).
Can someone please help me debug the issue?

Backtest Error

ValueError:Rank(A) < p or Rank([G; A]) < n

Above error happens when there are redundant constraints (and/or infinite solution choice) but I am unable to find any. Please help!

Below is the problem that I am trying to solve, it's math and code

Maximize c'w (maximize alpha)
where -0.01 < w < 0.01 (Position Concentration Limit)
0.925 < Sum(abs(w)) < 1.075 (Leverage/Gross Exposure limit)
-0.075 < Sum(w) < 0.075 (Net Exposure Limit)
Sum(abs(w-wo)) < 0.95 (Turnover Limit)
-0.375 < A'w < 0.375 (Style Exposure)
-0.175 < B'w < 0.175 (Sector Exposure)

Below, create new variables to convert all non-linear restrictions to linear restrictions

w = wL + wS (Rewriting portfolio weights as sum of long and short weights)
0 < wL < U
-U < wS < 0

Rewriting Leverage Constraint

lLeverage < Sum(wL) - Sum(wS) < uLeverage

Rewriting Net Exposure Constraint

lNetExposure < Sum(wL) + Sum(wS) < uNetExposure

Rewriting Turnover Constraint (creating new variable tw)

tw = abs(w - w0)
tw > w - w0
tw > - (w - w0)
tw > 0

Sum(tw) < MaxTurnover

Rewriting Style Exposure

Min Style Exposure < A'wL + A'wS < Max Style Exposure

Rewriting Sector Exposure

Min Sector Exposure < B'wL + B'wS < Max Sector Exposure

Python Code (More in attached backtest) 

def optimize(context):  
    alpha = context.alpha  
    factor_loadings = context.risk_factor_betas.loc[alpha.index.tolist(), :]  
    beta = context.beta_pipeline[['beta']].loc[alpha.index.tolist(), :]  
    nstocks = alpha.size  
    initial_portfolio = context.initial_portfolio  
    min_turnover = context.min_turnover  
    max_turnover = context.max_turnover - min_turnover  
    # w = wl + ws  
    # Number of variables =  
    #    2*nStocks(weights) +  
    #    nStocks(turnover link vars)  
    #Number of Inequality group restrictions =  
    #    2*nSector +  
    #    2*nStyles  +  
    #    2*1 (net exposure restrictions)  
    #    2*1 (gross exposure restriction) +  
    #    1 (turnover restriction)  
    #    2*nStocks(turnover link vars) +  
    #    2*1 (market beta exposure) +  
    #    2*nStocks (Weight long bounds)  
    #    2*nStocks (Weight short bounds)  
    #    nStocks (turnover link vars lower bound)  
    ones_matrix = matrix(1.0, (nstocks, 1), 'd')  
    zeros_matrix = matrix(0.0, (nstocks, 1), 'd')  

    # Group Constraints - 1/2  
    # min_exposure < Risk Loading transpose * W < max_exposure  
    sector_exp_matrix = matrix()  
    nsectors = len(context.sectors)  
    sector_exp_bounds = matrix(context.max_sector_exposure, (nsectors,1), 'd')

    for col_name in context.sectors:  
        _loadings = sparse([  
                matrix(factor_loadings[col_name].values, (nstocks, 1), 'd'),  
                matrix(factor_loadings[col_name].values, (nstocks, 1), 'd'),  
                zeros_matrix]).T  
        if sector_exp_matrix:  
            sector_exp_matrix = sparse([sector_exp_matrix, _loadings])  
        else:  
            sector_exp_matrix = _loadings  
    style_exp_matrix = matrix()  
    nstyles = len(context.styles)  
    style_exp_bounds = matrix(context.max_style_exposure, (nstyles,1), 'd')

    for col_name in context.styles:  
        _loadings = sparse([  
                matrix(factor_loadings[col_name].values, (nstocks, 1), 'd'),  
                matrix(factor_loadings[col_name].values, (nstocks, 1), 'd'),  
                zeros_matrix]).T  
        if style_exp_matrix:  
            style_exp_matrix = sparse([style_exp_matrix, _loadings])  
        else:  
            style_exp_matrix = _loadings  
    # Group Constraints - 3/4  
    # Dollar Neutral or min_net_exp < sum of weights < max_net_exp  
    net_exp_matrix = sparse([  
            ones_matrix,  
            ones_matrix,  
            zeros_matrix]).T  
    net_exp_bounds = matrix(context.max_net_exposure, (1,1), 'd')  
    # Group Constraints - 5/6  
    # Gross Exposure or min_gross_exp < sum of abs weights < max_gross_exp  
    gross_exp_matrix = sparse([  
            ones_matrix,  
            -1*ones_matrix,  
            zeros_matrix]).T  
    gross_exp_bounds_u = matrix(context.max_leverage, (1,1), 'd')  
    gross_exp_bounds_l = -1*matrix(context.min_leverage, (1,1), 'd')  
    # Group Constraints - 8  
    # Turnover < max_turnover  
    # Sum of abs(w-w0) < ma_turnover  
    # dw = abs(w-w0)  (new variable)  
    # dw >= w - wo  
    # dw >= -(w - wo)  
    # dw >= 0 variable constraint  
    # dw1 + dw2 + ... + dwn < max_turnover (Number of Group Restriction = 1)  
    # Total Number of Group Restriction = nFactors + nStocks  
    turnover_matrix = sparse([  
            zeros_matrix,  
            zeros_matrix,  
            ones_matrix]).T  
    turnover_bound = matrix(max_turnover, (1,1), 'd')  
    turnover_lnk_bounds = matrix(initial_portfolio, (nstocks,1), 'd')  
    # d = |wl + ws - w0|  
    # d > wl + ws - w0  => wl + ws - d < w0  
    # d > -wl -ws + w0 => -wl -ws -d < -w0  
    # Group Constraints - 9 (turnover link vars)  
    turnover_lnk_matrix_1 = matrix()  
    for i in range(nstocks):  
        var_indicator = np.zeros(nstocks)  
        var_indicator[i] = 1.0  
        _var_link_matrix =  sparse([  
                matrix(var_indicator, (nstocks,1), 'd'),  
                matrix(var_indicator, (nstocks,1), 'd'),  
                matrix(-var_indicator, (nstocks,1), 'd')]).T  
        if turnover_lnk_matrix_1:  
            turnover_lnk_matrix_1 = sparse([turnover_lnk_matrix_1, _var_link_matrix])  
        else:  
            turnover_lnk_matrix_1 = _var_link_matrix  

    # Group Constraints - 10 (turnover link vars)  
    turnover_lnk_matrix_2 = matrix()  
    for i in range(nstocks):  
        var_indicator = np.zeros(nstocks)  
        var_indicator[i] = 1.0  
        _var_link_matrix =  sparse([  
                matrix(-var_indicator, (nstocks,1), 'd'),  
                matrix(-var_indicator, (nstocks,1), 'd'),  
                matrix(-var_indicator, (nstocks,1), 'd')]).T  
        if turnover_lnk_matrix_2:  
            turnover_lnk_matrix_2 = sparse([turnover_lnk_matrix_2, _var_link_matrix])  
        else:  
            turnover_lnk_matrix_2 = _var_link_matrix  
    # Group constraints - 11/12  
    # Market beta exposure  
    # lBeta < (B1*W1 + B2*W2 + ... + BnWn) < uBeta  
    market_beta_exp_matrix = sparse([  
            matrix(beta.values, (nstocks,1), 'd'),  
            matrix(beta.values, (nstocks,1), 'd'),  
            zeros_matrix]).T  
    market_beta_exp_bound = matrix(context.max_beta, (1,1), 'd')  
    # Group Restrictions - 13/14  
    # Portfolio Weight Upper Bound  
    weight_var_matrix_long = matrix()  
    for i in range(nstocks):  
        var_indicator = np.zeros(nstocks)  
        var_indicator[i] = 1.0  
        _var_matrix =  sparse([  
                matrix(var_indicator, (nstocks,1), 'd'),  
                zeros_matrix,  
                zeros_matrix]).T  
        if weight_var_matrix_long:  
            weight_var_matrix_long = sparse([weight_var_matrix_long, _var_matrix])  
        else:  
            weight_var_matrix_long = _var_matrix  

    weight_var_matrix_short = matrix()  
    for i in range(nstocks):  
        var_indicator = np.zeros(nstocks)  
        var_indicator[i] = 1.0  
        _var_matrix =  sparse([  
                zeros_matrix,  
                matrix(var_indicator, (nstocks,1), 'd'),  
                zeros_matrix]).T  
        if weight_var_matrix_short:  
            weight_var_matrix_short = sparse([weight_var_matrix_short, _var_matrix])  
        else:  
            weight_var_matrix_short = _var_matrix  

    turnover_var_matrix = matrix()  
    for i in range(nstocks):  
        var_indicator = np.zeros(nstocks)  
        var_indicator[i] = -1.0  
        _var_matrix =  sparse([  
                zeros_matrix,  
                zeros_matrix,  
                matrix(var_indicator, (nstocks,1), 'd')]).T  
        if turnover_var_matrix:  
            turnover_var_matrix = sparse([turnover_var_matrix, _var_matrix])  
        else:  
            turnover_var_matrix = _var_matrix  
    # Combine all Group Restrictons  
    G = sparse([  
             sector_exp_matrix,  
             -1*sector_exp_matrix,  
               style_exp_matrix,  
             -1*style_exp_matrix,  
            net_exp_matrix,  
             -1*net_exp_matrix,  
             gross_exp_matrix,  
             -1*gross_exp_matrix,  
            turnover_matrix,  
            turnover_lnk_matrix_1,  
            turnover_lnk_matrix_2,  
            market_beta_exp_matrix,  
           -1*market_beta_exp_matrix,  
           -1*weight_var_matrix_long,  
            weight_var_matrix_long,  
            -1*weight_var_matrix_short,  
            weight_var_matrix_short,  
            turnover_var_matrix  
        ])  


    h = matrix([  
            sector_exp_bounds,  
            sector_exp_bounds,  
            style_exp_bounds,  
            style_exp_bounds,  
            net_exp_bounds,  
            net_exp_bounds,  
            gross_exp_bounds_u,  
            gross_exp_bounds_l,

            turnover_bound,  
            turnover_lnk_bounds,  
            -1*turnover_lnk_bounds,  
            market_beta_exp_bound,  
            market_beta_exp_bound,

            zeros_matrix, # xL > 0  
            matrix(context.max_pos_size, (nstocks, 1), 'd'), #xL < U  
            matrix(context.max_pos_size, (nstocks, 1), 'd'), #-xS < U  
            zeros_matrix, #xS < 0  
            zeros_matrix  
        ])  
    #Objective Function - Maximize Alpha  
    c = matrix([  
        matrix(alpha.values, (nstocks,1)),  
        matrix(alpha.values, (nstocks,1)),  
        zeros_matrix])  
    sol = solvers.lp(-c, G, h)  
    print sol['x']
3 responses
Shiv Chawla
May 23, 2019

Surprisingly, it works for nstocks = 500 but fails after that. This is strange!

Shiv Chawla
May 24, 2019

Thanks Satya for the help! CVXPY makes optimization very easy. Brilliant!!

I have modeled the problem using CVXPY and could solve the problem (BUT turnover restrictions are throwing a curveball). With absolute restrictions turned to two linear restrictions, the optimization setup cannot be solved by CVXPY (using ECOS). Could you please take a look?

Turnover Restrictions

t = |wL + wS - w0|
t > wL + wS - w0 => -np.inf < wL + wS - t < w0
t > -wL -wS + w0 => w0 < wL + wS + t < np.inf
t >= 0

Optimization Code (more in backtest)

def optimize(context, alpha, factor_loadings, beta):  
    nstocks = alpha.size  
    initial_portfolio = context.initial_portfolio.reshape((nstocks,))  
    min_turnover = context.min_turnover  
    max_turnover = context.max_turnover - min_turnover  
    # w = wl + ws  
    # Number of variables =  
    #    2*nStocks(weights) +  
    #    nStocks(turnover link vars)  
    #Number of Inequality group restrictions =  
    #    nSector +  
    #    nStyles  +  
    #    1 (net exposure restrictions)  
    #    1 (gross exposure restriction) +  
    #    1 (turnover restriction)  
    #    2*nStocks(turnover link vars) +  
    #    1 (market beta exposure) +  

    zeros_array = np.zeros(nstocks)  
    I = np.identity(nstocks)  
    ONE_ROW_MAT = np.asmatrix(np.full((1, nstocks), 1.0))  
    ZERO_ROW_MAT = np.asmatrix(np.full((1, nstocks), 0.0))  
    # Group Constraints - 1  
    # min_exposure < Risk Loading transpose * W < ma_exposure  
    sector_exp_matrix = None  
    nsectors = len(context.sectors)  
    sector_exp_bounds = np.full((nsectors), context.max_sector_exposure)  

    for col_name in context.sectors:  
        _loadings = np.matrix(np.hstack((  
                factor_loadings[col_name].values,  
                factor_loadings[col_name].values,  
                zeros_array)).reshape((1, 3*nstocks)))  
        if sector_exp_matrix is not None:  
            sector_exp_matrix = np.concatenate((sector_exp_matrix, _loadings))  
        else:  
            sector_exp_matrix = _loadings  
    sector_exp_matrix = scipysp.csc_matrix(sector_exp_matrix)  
    style_exp_matrix = None  
    nstyles = len(context.styles)  
    style_exp_bounds = np.full((nstyles), context.max_style_exposure)

    for col_name in context.styles:  
        _loadings = np.matrix(np.hstack((  
                factor_loadings[col_name].values,  
                factor_loadings[col_name].values,  
                zeros_array)).reshape((1, 3*nstocks)))  
        if style_exp_matrix is not None:  
            style_exp_matrix = np.concatenate((style_exp_matrix, _loadings))  
        else:  
            style_exp_matrix = _loadings  
    style_exp_matrix = scipysp.csc_matrix(style_exp_matrix)  
    # Group Constraints - 2  
    # Dollar Neutral or min_net_exp < sum of weights < max_net_exp  
    net_exp_matrix = scipysp.csc_matrix(  
        np.concatenate((ONE_ROW_MAT, ONE_ROW_MAT, ZERO_ROW_MAT), axis=1))  
    net_exp_bounds = np.full((1), context.max_net_exposure)  
    # Group Constraints - 3  
    # Gross Exposure or min_gross_exp < sum of abs weights < max_gross_exp  
    # aw = abs(w)  
    # w - aw <= 0 ; aw > 0  
    # min_gross_exp < aw1 + aw2 + ..... + awn < max_gross_exp  
    gross_exp_matrix = scipysp.csc_matrix(  
        np.concatenate((ONE_ROW_MAT, -ONE_ROW_MAT, ZERO_ROW_MAT), axis=1))  
    gross_exp_bounds_u = np.full((1), context.max_leverage)  
    gross_exp_bounds_l = np.full((1), context.min_leverage)  
    # Group Constraints - 4  
    # Turnover < max_turnover  
    # Sum of abs(w-w0) < ma_turnover  
    # t = abs(w-w0) #new variable  
    turnover_matrix = scipysp.csc_matrix(  
        np.concatenate((ZERO_ROW_MAT, ZERO_ROW_MAT, ONE_ROW_MAT), axis=1))  
    turnover_bound = np.full((1), max_turnover)  
    # t = |wl + ws - w0|  
    # t > wl + ws - w0  =>  -np.inf < wl + ws - t < w0  
    # t > -wl -ws + w0 =>  w0 < wl + ws + t < np.inf  
    # t >= 0  
    # Group Constraints - 5 (turnover link vars)  
    turnover_lnk_matrix_1 = scipysp.csc_matrix(  
        np.concatenate((I, I, -1*I), axis=1))  
    turnover_lnk_matrix_2 = scipysp.csc_matrix(  
        np.concatenate((I, I, I), axis=1))  
    turnover_lnk_bounds_1l = -np.full((nstocks), np.inf)  
    turnover_lnk_bounds_1u = initial_portfolio  
    turnover_lnk_bounds_2l = initial_portfolio  
    turnover_lnk_bounds_2u = np.full((nstocks), np.inf)

    # Group constraints - 6  
    # Market beta exposure  
    # lBeta < (B1*W1 + B2*W2 + ... + BnWn) < uBeta  
    market_beta_exp_matrix = scipysp.csc_matrix(np.matrix(np.concatenate((  
            np.reshape(beta.values, (nstocks,)),  
            np.reshape(beta.values, (nstocks,)),  
            zeros_array)).reshape((1, 3*nstocks))))  
    market_beta_exp_bound = np.full((1), context.max_beta)

    # Create optimization variables  
    wL = cp.Variable(nstocks)  
    wS = cp.Variable(nstocks)  
    tr = cp.Variable(nstocks)

    # Combine all Group Restrictons  
    A = cp.vstack(  
            sector_exp_matrix,  
            style_exp_matrix,  
            net_exp_matrix,  
            gross_exp_matrix,  
            turnover_matrix,  
            # turnover_lnk_matrix_1,  
            # turnover_lnk_matrix_2,  
            market_beta_exp_matrix,  
        )  
    # Group Restrictions Upper Bound  
    Ub = np.hstack((  
            sector_exp_bounds,  
            style_exp_bounds,  
            net_exp_bounds,  
            gross_exp_bounds_u,  
            turnover_bound,  
            # turnover_lnk_bounds_1u,  
            # turnover_lnk_bounds_2u,  
            market_beta_exp_bound  
        ))  
    # Group Restrictions Lower Bound  
    Lb = np.hstack((  
            -1*sector_exp_bounds,  
            -1*style_exp_bounds,  
            -1*net_exp_bounds,  
            gross_exp_bounds_l,  
            np.full((1,), 0.0),  
            # turnover_lnk_bounds_1l,  
            # turnover_lnk_bounds_2l,  
            -1*market_beta_exp_bound  
        ))  
    # Stack Variables  
    x = cp.vstack(wL,wS,tr)  
    # Optimization Problem Constraints (Group + Variable)  
    constraints = [  
        A*x <= Ub,  
        A*x >= Lb,  
        wL >= 0,  
        wL <= context.max_pos_size,  
        wS <= 0,  
        wS >= -context.max_pos_size,  
        tr >= 0]  
    #Objective Function - Maximize Alpha  
    c = np.hstack((  
        alpha.values,  
        alpha.values,  
        zeros_array))  
    objective = cp.Maximize(c.T*x)  
    prob = cp.Problem(objective, constraints)  
    prob.solve()
Shiv Chawla
May 24, 2019

Fixed the turnover constraint but no success. Still the solver issue!!
Fixed Turnover Constraint

t = |wL + wS - w0|
t1 > wL + wS - w0 => -np.inf < wL + wS - t1 < w0
t2 > -wL -wS + w0 => w0 < wL + wS + t2 < np.inf
t1 >= 0; t2 >=0

But then I used all the features of CVXPY like norm and sum_entries to replace net/gross/turnover restrictions and finally some success!!
However, the optimized portfolio doesn't match the Q optimized portfolio (I guess the style/sector restrictions are different).

Are the style/sector limits optimization limits mentioned anywhere in the documentation?

Also, Does Q investment team use any minimum trade size limits? Optimization may generate some very small trades which may be not be desirable.

About Quantopian
Careers
Community
Events
Documentation
Academia
Lectures
Status
Facebook
Twitter
LinkedIn
YouTube
Terms of Use
Privacy Policy

The material on this website is provided for informational purposes only and does not constitute an offer to sell, a solicitation to buy, or a recommendation or endorsement for any security or strategy, nor does it constitute an offer to provide investment advisory services by Quantopian.

In addition, the material offers no opinion with respect to the suitability of any security or specific investment. No information contained herein should be regarded as a suggestion to engage in or refrain from any investment-related course of action as none of Quantopian nor any of its affiliates is undertaking to provide investment advice, act as an adviser to any plan or entity subject to the Employee Retirement Income Security Act of 1974, as amended, individual retirement account or individual retirement annuity, or give advice in a fiduciary capacity with respect to the materials presented herein. If you are an individual retirement or other investor, contact your financial advisor or other fiduciary unrelated to Quantopian about whether any given investment idea, strategy, product or service described herein may be appropriate for your circumstances. All investments involve risk, including loss of principal. Quantopian makes no guarantees as to the accuracy or completeness of the views expressed in the website. The views are subject to change, and may have become unreliable for various reasons, including changes in market conditions or economic circumstances.