Put/Call efficient frontier model#

PutCall.py WorldIndices.gdx

"""
## GAMSSOURCE: https://www.gams.com/latest/finlib_ml/libhtml/finlib_PutCall.html
## LICENSETYPE: Demo
## MODELTYPE: LP
## DATAFILES: WorldIndices.gdx


Put/Call efficient frontier model

MAD.gms: Put/Call efficient frontier model.
Consiglio, Nielsen and Zenios.
PRACTICAL FINANCIAL OPTIMIZATION: A Library of GAMS Models, Section 5.7
Last modified: Apr 2008.
"""

from __future__ import annotations

import os
from pathlib import Path
from sys import argv

import numpy as np
import pandas as pd
from gamspy import (
    Alias,
    Card,
    Container,
    Equation,
    Model,
    Number,
    Options,
    Parameter,
    Set,
    Sum,
    Variable,
)


def index_data():
    data = np.array(
        [
            1.034769211,
            1.024362083,
            1.005721076,
            1.014359531,
            1.008024402,
            1.015164086,
            1.011776948,
            0.999367312,
            1.022337255,
            1.014084294,
            1.002478561,
            1.004830772,
            1.022048865,
            1.001584376,
            1.014789245,
            1.018465908,
            1.02023507,
            1.001142139,
            1.028357519,
            1.012274672,
            1.005320247,
            0.991451365,
            1.007344745,
            1.017165464,
            0.997931432,
            0.990041853,
            0.9933217,
            0.991912481,
            1.036980768,
            0.994451813,
            1.007412617,
            0.955150659,
            0.958840232,
            1.033624988,
            1.018079311,
            1.010736753,
            1.016700717,
            1.048121923,
            1.015146812,
            1.0074283,
            1.025778693,
            0.983151788,
            1.022476292,
            1.004850227,
            1.0091122,
            0.995361207,
            0.992432266,
            1.024212556,
            1.022181138,
            1.006772127,
            0.990354861,
            1.0067649,
            1.008943286,
            0.983974563,
            1.000347655,
            1.001230771,
            1.020880921,
            1.019415483,
            1.008638044,
            1.017173042,
            1.011831605,
            1.018260066,
            1.014054495,
            1.011841085,
            0.988155751,
            1.009224769,
            1.010321982,
            1.033008126,
            0.999626315,
            1.022733809,
            1.00083866,
            1.039925522,
            1.0144168,
            0.974402144,
            0.971108237,
            1.002198078,
            1.001221179,
            0.988342113,
            1.01750901,
            1.019968377,
            0.992895551,
            1.000905761,
            0.997333999,
            0.994982525,
            0.987972315,
            1.000443735,
            1.01583315,
            1.018338346,
            1.026183302,
            0.999209138,
            1.018809687,
            1.009775599,
            1.01722367,
            0.999752239,
            1.015986198,
            1.019400469,
            1.017624299,
            0.996067065,
            1.00947457,
            1.0164916,
            0.997929107,
            1.008892447,
            0.990498799,
            1.011843222,
            1.023338977,
            1.001536085,
            1.023978065,
            0.994238826,
            1.022819962,
            1.012142345,
            0.989862601,
            1.021050206,
            1.019165002,
            1.023260624,
            1.027413626,
            0.971346297,
            1.026087902,
            0.959990261,
            1.008484984,
            1.014365951,
            1.010347985,
            1.029266037,
            1.01989954,
            0.999383722,
            0.992191473,
            0.993782093,
            1.010730888,
            0.918149916,
            1.008729179,
            1.031641438,
            1.032632304,
            1.004115995,
            1.007398109,
            0.99859111,
            1.032604061,
            1.027883261,
            0.979268923,
            1.024656356,
            0.989556024,
            1.001926578,
            0.991713793,
            1.020697126,
            1.022268183,
            1.031273723,
            0.992301368,
            1.011547524,
            1.02121809,
            0.988701862,
            0.991540642,
            1.014546705,
            0.997009882,
            1.01168258,
            0.982543898,
            0.992197467,
        ]
    )
    return data


def main():
    gdx_file = str(Path(__file__).parent.absolute()) + "/WorldIndices.gdx"
    m = Container(
        system_directory=os.getenv("SYSTEM_DIRECTORY", None),
        load_from=gdx_file,
    )

    output = argv[1] if len(argv) > 1 else "PutCallModel"

    # SETS #
    i, l = m.getSymbols(["i", "l"])

    # SCALARS #
    Budget = Parameter(
        m, name="Budget", description="Nominal investment budget"
    )
    Omega = Parameter(
        m, name="Omega", description="Bound on the expected shortfalls"
    )

    Budget[...] = 100.0

    # PARAMETERS #
    pr = Parameter(m, name="pr", domain=l, description="Scenario probability")
    P = Parameter(m, name="P", domain=[i, l], description="Final values")
    EP = Parameter(m, name="EP", domain=i, description="Expected final values")
    TargetIndex = Parameter(
        m, name="TargetIndex", domain=l, description="Target index returns"
    )
    AssetReturns = m.getSymbols(["AssetReturns"])[0]

    pr[l] = 1.0 / Card(l)

    P[i, l] = 1 + AssetReturns[i, l]

    EP[i] = Sum(l, pr[l] * P[i, l])

    # To test the model with a market index, uncomment the following line two lines.
    # Note that, this index is consistent only when using WorldIndexes.inc.
    Index = Parameter(
        m,
        name="Index",
        domain=l,
        records=index_data(),
        description="Index returns",
    )
    TargetIndex[l] = Index[l]

    # VARIABLES
    yPos = Variable(
        m,
        name="yPos",
        type="positive",
        domain=l,
        description="Positive deviations",
    )
    yNeg = Variable(
        m,
        name="yNeg",
        type="positive",
        domain=l,
        description="Negative deviations",
    )
    x = Variable(
        m,
        name="x",
        type="free",
        domain=i,
        description="Holdings of assets in monetary units (not proportions)",
    )

    # EQUATIONS #
    BudgetCon = Equation(
        m,
        name="BudgetCon",
        type="regular",
        description="Equation defining the budget contraint",
    )
    TargetDevDef = Equation(
        m,
        name="TargetDevDef",
        type="regular",
        domain=l,
        description="Equations defining the positive and negative deviations",
    )
    PutCon = Equation(
        m,
        name="PutCon",
        type="regular",
        description=(
            "Constraint to bound the expected value of the negative deviations"
        ),
    )

    BudgetCon[...] = Sum(i, x[i]) == Budget

    PutCon[...] = Sum(l, pr[l] * yNeg[l]) <= Omega

    TargetDevDef[l] = (
        Sum(i, (P[i, l] - TargetIndex[l]) * x[i]) == yPos[l] - yNeg[l]
    )

    # Objective function definition for MAD
    ObjDef = Sum(l, pr[l] * yPos[l])

    UnConPutCallModel = Model(
        m,
        name="UnConPutCallModel",
        equations=[PutCon, TargetDevDef],
        problem="LP",
        sense="MAX",
        objective=ObjDef,
    )

    # Set the average level of downside (risk) allowed

    Omega[...] = 0.1
    UnConPutCallModel.solve(
        options=Options(
            variable_listing_limit=0,
            equation_listing_limit=0,
            solver_link_type=2,
        )
    )

    # Dual of the UnConstrained Put/Call model

    # VARIABLES #
    Pi = Variable(m, name="Pi", type="positive", domain=l)
    PiOmega = Variable(m, name="PiOmega", type="positive")

    # EQUATIONS #
    DualTrackingDef = Equation(
        m, name="DualTrackingDef", type="regular", domain=i
    )
    MeasureDef = Equation(m, name="MeasureDef", type="regular", domain=l)

    DualTrackingDef[i] = Sum(l, (P[i, l] - TargetIndex[l]) * Pi[l]) == 0.0

    MeasureDef[l] = pr[l] * PiOmega - Pi[l] >= 0

    DualObjDef = Omega * PiOmega  # Objective function definition Dual

    Pi.lo[l] = pr[l]

    UnConDualPutCallModel = Model(
        m,
        name="UnConDualPutCallModel",
        equations=[DualTrackingDef, MeasureDef],
        problem="LP",
        sense="MIN",
        objective=DualObjDef,
    )

    UnConDualPutCallModel.solve(
        options=Options(
            variable_listing_limit=0,
            equation_listing_limit=0,
            solver_link_type=2,
        )
    )

    # Display PiOmega.l and Pi.l and check that they are, respectively, equal
    # to TargetDevDef.m and PutCon.m
    # GAMS provides the dual prices directly, so it is not
    # really necessary to build explicitly the dual model.

    # PARAMETER
    PrimalDual = Parameter(
        m,
        name="PrimalDual",
        domain=[l, "*"],
        description="Compare primal and dual soultions",
    )

    PrimalDual[l, "pi.l"] = -Pi.l[l]
    PrimalDual[l, "TargetDevDef.m"] = TargetDevDef.m[l]
    PrimalDual[l, "Difference"] = TargetDevDef.m[l] + Pi.l[l]

    # DISPLAY UnConDualPutCallModel.objective_value,PiOmega.l,PutCon.m,PrimalDual

    # We propose an alternative way to build a frontier using
    # the loop statement. Such a structure is suitable for the
    # GDX utility (for details, see gdxutility.pdf included in the doc folder )

    # SET
    FrontierPoints = Set(
        m,
        name="FrontierPoints",
        records=[f"P_{fp}" for fp in range(1, 51)],
        description="Number of points in the frontier",
    )
    j = Alias(m, name="j", alias_with=FrontierPoints)

    # PARAMETER
    FrontierPortfolios = Parameter(
        m,
        name="FrontierPortfolios",
        domain=[j, i],
        description="Frontier portfolios",
    )
    CallValues = Parameter(
        m, name="CallValues", domain=[j, "*"], description="Call values"
    )
    DualPrices = Parameter(
        m, name="DualPrices", domain=[j, "*"], description="Dual prices"
    )
    PutCall = Parameter(
        m, name="PutCall", domain=[j, "*"], description="Put and Call values"
    )
    OmegaLevels = Parameter(
        m, name="OmegaLevels", domain=j, description="Risk levels (Omega)"
    )

    # We assign to each point a risk level Omega

    OmegaLevels["P_1"] = 0.01

    for idx, jj in enumerate(j.toList()):
        if idx == 0:
            continue
        elif (idx + 1) <= 10:
            OmegaLevels[jj] = OmegaLevels[j.toList()[idx - 1]] + Number(0.01)
        elif (idx + 1) > 10:
            OmegaLevels[jj] = OmegaLevels[j.toList()[idx - 1]] + Number(0.025)

    DualPrices[j, "Omega"] = OmegaLevels[j]

    CallValues[j, "Omega"] = OmegaLevels[j]

    # Set some liquidity constraints

    x.lo[i] = -100.0
    x.up[i] = 100.0

    for jj in j.toList():
        Omega[...] = OmegaLevels[jj]
        UnConPutCallModel.solve(
            options=Options(
                variable_listing_limit=0,
                equation_listing_limit=0,
                solver_link_type=2,
            )
        )

        FrontierPortfolios[jj, i] = x.l[i]
        CallValues[jj, "Mild Constraint"] = UnConPutCallModel.objective_value
        DualPrices[jj, "Mild Constraint"] = PutCon.m

    # EXPORT RESULT SUMMARY TO EXCEL SHEET #
    writer = pd.ExcelWriter(f"{output}.xlsx", engine="openpyxl")
    FrontierPortfolios.pivot().to_excel(writer, sheet_name="MildPortfolios")

    # Set tight liquidity constraints

    x.lo[i] = -20.0
    x.up[i] = 20.0

    for jj in j.toList():
        Omega[...] = OmegaLevels[jj]
        UnConPutCallModel.solve(
            options=Options(
                variable_listing_limit=0,
                equation_listing_limit=0,
                solver_link_type=2,
            )
        )

        FrontierPortfolios[jj, i] = x.l[i]
        CallValues[jj, "Tight Constraint"] = UnConPutCallModel.objective_value
        DualPrices[jj, "Tight Constraint"] = PutCon.m

    DualPrices.pivot().to_excel(writer, sheet_name="DualPrices")
    CallValues.pivot().to_excel(writer, sheet_name="CallValues")
    FrontierPortfolios.pivot().to_excel(writer, sheet_name="TightPortfolios")

    # Determine the liquidity and discount premium
    # for one put/call efficient portfolio.

    Omega[...] = 0.475

    UnConPutCallModel.solve(
        options=Options(
            variable_listing_limit=0,
            equation_listing_limit=0,
            solver_link_type=2,
        )
    )

    # SCALAR
    Df = Parameter(m, name="Df")

    # PARAMETERS
    Price = Parameter(m, name="Price", domain=i)
    Discount = Parameter(m, name="Discount", domain=i)
    Premium = Parameter(m, name="Premium", domain=i)
    BenchMarkNeutralPrice = Parameter(
        m, name="BenchMarkNeutralPrice", domain=i
    )
    Psi = Parameter(m, name="Psi", domain=l)
    liquidity = Parameter(
        m, name="liquidity", domain=[i, "*"], description="Liquidity report"
    )

    Discount[i] = 0.0

    Premium[i] = 0.0

    Df[...] = Sum(l, -TargetDevDef.m[l])

    Psi[l] = -TargetDevDef.m[l] / Df

    BenchMarkNeutralPrice[i] = Sum(l, Psi[l] * P[i, l]) / Sum(
        l, Psi[l] * TargetIndex[l]
    )

    Discount[i].where[x.m[i] > 0] = x.m[i] / Sum(
        l, (-TargetDevDef.m[l]) * TargetIndex[l]
    )

    Premium[i].where[x.m[i] < 0] = (-x.m[i]) / Sum(
        l, (-TargetDevDef.m[l]) * TargetIndex[l]
    )

    Price[i] = BenchMarkNeutralPrice[i] + Premium[i] - Discount[i]

    liquidity[i, "Premium"] = Premium[i]
    liquidity[i, "Discount"] = Discount[i]

    liquidity.pivot().to_excel(writer, sheet_name="Liquidity")

    # Put/call model with balance constraint

    PutCallModel = Model(
        m,
        name="PutCallModel",
        equations=[BudgetCon, PutCon, TargetDevDef],
        problem="LP",
        sense="MAX",
        objective=ObjDef,
    )

    for jj in j.toList():
        Omega[...] = OmegaLevels[jj]
        PutCallModel.solve(
            options=Options(
                variable_listing_limit=0,
                equation_listing_limit=0,
                solver_link_type=2,
            )
        )
        FrontierPortfolios[jj, i] = x.l[i]
        PutCall[jj, "Put side"] = PutCon.l
        PutCall[jj, "Call side"] = PutCallModel.objective_value

    PutCall.pivot().to_excel(writer, sheet_name="Frontiers")
    FrontierPortfolios.pivot().to_excel(writer, sheet_name="Portfolios")
    writer.close()


if __name__ == "__main__":
    main()