# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: Python 3
#     language: python
#     name: python3
# ---

# # ESG-score using a MiniMax Model

# ## Problem setting:

# As a Business-idea firm xy plans to launch a new app, in which clients are able to invest in a socially responsible manner. The app has to be able to compute for every portfolio an Environmental, Social and Governance (ESG) score between 0 and 100, where a higher score indicates a more responsible investment. Using an optimization process for each client a customized portfolio that takes into account the provided target values should be derived. Since firm xy is based in Switzerland, the initial investments universe comprises the stocks listed on SMI. The files prices1617.csv, prices18.csv and scores.csv contain the weekly closing prices of the SMI and its constituents for the years 2016-2017 and 2018 as well as the latest ESG scores for all SMI constituents, respectively.

# Import packages
import pandas as pd
import numpy as np
import gurobipy as gb
import matplotlib.pyplot as plt

# Define precision of number display
pd.set_option('precision', 5)

# Import prices
weekly_prices = pd.read_csv('prices1617.csv', index_col='date')
# Look at first five rows
weekly_prices.head()

# Compute weekly returns
weekly_returns = weekly_prices / weekly_prices.shift(1) - 1
# Drop NA-return (first period)
weekly_returns = weekly_returns.dropna()
# Look at first five rows
weekly_returns.head()

# Compute expected returns
all_expected_returns = weekly_returns.mean()
all_expected_returns

# Get array of stocks (without index SMI)
stocks = weekly_returns.columns.drop('SMI')

# ## Portfolio 1 (P1):
# As a first idea for a portfolio selection an equally weighted portfolio approach is chosen. This portfolio is referred by P1.
# In a first step the variance, mean-absolute deviation, Beta (relative to SMI), minimum return and the expected downside $(r_t =1\%)$ for the portfolio P1 as well as for the SMI are computed.

# Create equally weighted portfolio
weights_p1 = pd.Series(1 / len(stocks), index=stocks)
weekly_returns['P1'] = weekly_returns[stocks].dot(weights_p1)

# Calculate variance
weekly_returns[['P1', 'SMI']].var()

# Calculate mean-absolute deviation
weekly_returns[['P1', 'SMI']].mad()

# Compute the beta of P1
weekly_returns['P1'].cov(weekly_returns['SMI']) / weekly_returns['SMI'].var()

# Compute the beta of SMI
weekly_returns['SMI'].cov(weekly_returns['SMI']) / weekly_returns['SMI'].var()

# Compute the minimum return
weekly_returns[['P1', 'SMI']].min()

# +
# Define target return
target_return = 0.01

# Compute expected downside with respect to given target return
np.maximum(0, -(weekly_returns[['P1', 'SMI']] - target_return)
          ).sum() / weekly_returns['P1'].size
# -

# ### Note:
# In a next step the ESG Environmental, ESG Governance and ESG Social scores of P1 are computed.

# Import ratings
scores = pd.read_csv('scores.csv', index_col='scores')

# Look at first five rows
scores.head()

# Get ESG scores for P1
scores.dot(weights_p1)

# ### Question:
# Does a distinct relationship between the expected returns and the average ESG scores exist? Do more sustainable stocks tend to have greatly inferior performance than less sustainable stocks?

# Get average ESG scores
average_esg_score = scores.mean()

# Visualize average ESG scores and expected returns in scatterplot
plt.figure(figsize=(9, 3), dpi=120)
plt.scatter(average_esg_score, all_expected_returns[stocks], alpha=0.6)
plt.xlabel('ESG average score')
plt.ylabel('Expected return');

# The data does not support the claim. It appears that no distinct relationship between the expected returns and the average ESG scores exist.

# ## Portfolio 2 (P2):
# As an alternative the firm xy has the idea to implement a MiniMax model as the underlying portfolio selection model (P2). Further a blacklist is introduced that consist of two pharmaceutical companies: Novartis and Roche. The weights for stocks on the blacklist must be zero in all subsequent portfolios.
#

# ### Note:
# A MiniMax(MM) model trades off the minimum portfolio rerun against the expected portfolio return. Therefore the minimum portfolio return is a downside risk measure. It is possible to formulate a MM model as a linear optimization problem. 
#
# Some technical characteristics of the MM model:
# * Minimum return is homogeneous
# * Minimum return is monotone
# * Minimum return is subadditive
# * Minimum return is translational invariant
#
# => Therefore the Minimum return is coherent
#
#
#
# Notation: 
# * $\alpha$ Minimum portfolio return over all scenarios
#
# \begin{align}
# Max. \alpha\\ 
# \sum_{i=n}^n r_{il}x_i \geq \alpha \quad (l \in \Omega)\\ 
# \sum_{i=n}^n \bar r_ix_i = \mu \\
# \sum_{i=n}^n x_i =1 \\
# x_i \geq 0 \quad (i =1,...,n) \\
# \alpha \in \Re
# \end{align}

# Get array of scenarios
scenarios = weekly_returns.index
# Compute expected returns
expected_returns = all_expected_returns[stocks]
# Define target portfolio return
target_return = all_expected_returns['SMI']

# +
# Create model
m = gb.Model('MM')

# Add a variable for each stock
x = pd.Series(m.addVars(stocks), index=stocks)

# Add a variable which represents the minimum portfolio return
minimum_return = m.addVar(lb=weekly_returns[stocks].min().min())

# Define portfolio return
portfolio_return = expected_returns.dot(x)

# Add objective function to model
m.setObjective(minimum_return, gb.GRB.MAXIMIZE)

# Add minimum return constraints
for l in scenarios:
    m.addConstr(weekly_returns.loc[l, stocks].dot(x) >= minimum_return)
    
# Add budget constraint
m.addConstr(x.sum() == 1);

# Add target return constraint
target_return_constr = m.addConstr(portfolio_return == target_return)

# Add blacklist constraints
m.addConstr(x.loc['ROCHE'] == 0)
m.addConstr(x.loc['NOVARTIS'] == 0);

# Specify Gurobi options
m.setParam('OutputFlag', 0)

# Run optimization
m.optimize()

# Get weights of stocks in optimal solution
weights_p2 = pd.Series([var.X for var in x], index=stocks)

# Get ESG scores for P2
scores[stocks].dot(weights_p2)
# -

# ### Question:
# What is the minimum return of P2 and how does it stand up compared to the SMI?

# Get minimum return of P2
minimum_return_p2 = minimum_return.X
minimum_return_p2

# Is minimum return of P2 bigger than minimum return of SMI?
minimum_return.X > weekly_returns['SMI'].min()

# The minimum return of P2 (-0.03953) is larger (less negative) than the minimum return of the SMI (-0.06360) and thus, the SMI is not efficient

# Get weights of stocks in optimal solution
weights_p2

# ## Note:
# Adjust the target return constraint of the model such that the optimal portfolio has an expected return equal to or greater than the target return (P3). 
#
# \begin{align}
# Max. \alpha\\ 
# \sum_{i=n}^n r_{il}x_i \geq \alpha \quad (l \in \Omega)\\ 
# \sum_{i=n}^n \bar r_ix_i \geq \mu \\
# \sum_{i=n}^n x_i =1 \\
# x_i \geq 0 \quad (i =1,...,n) \\
# \alpha \in \Re
# \end{align}

# Remove old target return constraint
m.remove(target_return_constr)

# Add new target return constraint
target_return_constr = m.addConstr(portfolio_return >= target_return)

# Run optimization
m.optimize()

# Get minimum return of P3
minimum_return.X

# The minimum return impoves from -3.95% in P2 to -3.2% in P3.

# Get weights of stocks in optimal solution
weights_p3 = pd.Series([var.X for var in x], index=stocks)
weights_p3

# ## Portfolio 4 (P4)
# In order to develop a variant of the minimax model that additionally takes into account the investors' target values for the three ESG score a separate penalty for each of the three ESG scores is implemented. If the portfolio score is equal to or greater than the target value, then the penalty is zero. If the portfolio score is smaller than the target value, then the penalty corresponds to the difference between the target value and the portfolio score.
# $penalty_{ESG} = max(taregt ESG score - portfolio ESG score,0)$ 
#
# However, firm xy decides that a linear increasing penalty is not appropriate. It should be ensured that large shortfalls are penalized heavily. The structure of the three penalty functions can be seen in the following graphs:
#
# ![title](picture.PNG)
#
#
# Additionally the already introduced constraints have to hold:
# * The portfolio return does not fall below a prespecified minimum portfolio return $\alpha$ in any scenario.
# * The expected portfolio return is equal to or greater than a prespecified target return µ.
# * The sum of the weights is equal to 1.
# * All portfolio weights are non-negative.
# * The weights for stocks on the blacklist must be zero.
#
# Formulate a linear optimization problem in which the sum of the three penalties (given the structures in the graphs) is minimized subject to all constraints.

# ### Note:
# The resulting optimization problem:
#
# Notation:
# * $S_i^E$: ESG Envionmental score of stock i
# * $S_i^S$: ESG Social score of stock i
# * $S_i^G$: ESG Governance score of stock i
# * $S_t^E$: Target ESG Environmental value
# * $S_t^S$: Target ESG Socail value
# * $S_t^G$: Target ESG Governance value
# * $x_i$: weight of asset i
# * $\delta_{Ej}$:  Shortfall of portfolio ESG Environmental score with respect to target ESG Environmental value in jth segment of penalty function (j = 1; 2)
# * $\delta_{Sj}$: Shortfall of portfolio ESG Social score with respect to target ESG Social value in jth segment of penalty function (j = 1; 2)
# * $\delta_{Gj}$: Shortfall of portfolio ESG Governance score with respect to target ESG Governance value in jth segment of penalty function (j = 1; 2) 
#
#
# \begin{align}
# Min. (8/10)\delta_{E1} + (12/10)\delta_{E2} + (10/10)\delta_{S1} + (20/10)\delta_{S2} + (4/10) \delta_{G1} + (8/10)\delta_{G2}\\
# \sum_{i=n}^n r_{il}x_i \geq \alpha \quad (l \in \Omega)\\ 
# \sum_{i=n}^n \bar r_ix_i \geq \mu \\
# \sum_{i=n}^n x_i =1 \\
# \sum_{i=n}^n S_i^E x_i  \geq S_t^E - \delta_{E1} - \delta_{E2}\\
# \sum_{i=n}^n S_i^S x_i  \geq S_t^S - \delta_{S1} - \delta_{S2}\\
# \sum_{i=n}^n S_i^G x_i  \geq S_t^G - \delta_{G1} - \delta_{G2}\\
# \delta_{E1} \leq 10 \\
# \delta_{S1} \leq 10 \\
# \delta_{G1} \leq 10 \\
# x_i \geq 0 \quad (i =1,...,n) \\
# x_i = 0 \quad (i = Novartis, Roche)\\
# \delta_{E1},\delta_{E2},\delta_{S1},\delta_{S2},\delta_{G1},\delta_{G2} \geq 0\\
# \end{align}
#

# Define target ESG values
S_t_E = 100
S_t_S = 100
S_t_G = 100

# +
# Create model
m = gb.Model('MM2')

# Add a variable for each stock
x = pd.Series(m.addVars(stocks), index=stocks)

# Add variables for ESG environmental penalties
delta_e1 = m.addVar()
delta_e2 = m.addVar()

# Add variables for ESG social penalties
delta_s1 = m.addVar()
delta_s2 = m.addVar()

# Add variables for ESG governance penalties
delta_g1 = m.addVar()
delta_g2 = m.addVar()

# Define minimum return
minimum_return = minimum_return_p2

# Define portfolio return
portfolio_return = expected_returns.dot(x)

# Define portfolio environmental score
portfolio_environmental_score = scores.loc['ESG Environmental', stocks].dot(x)

# Define portfolio social score
portfolio_social_score = scores.loc['ESG Social', stocks].dot(x)

# Define portfolio governance score
portfolio_governance_score = scores.loc['ESG Governance', stocks].dot(x)

# Add objective function to model
m.setObjective((8 / 10) * delta_e1 + (12 / 10) * delta_e2
               + (10 / 10) * delta_s1 + (20 / 10) * delta_s2
               + (4 / 10) * delta_g1 + (8 / 10) * delta_g2, gb.GRB.MINIMIZE)

# Add constraints for ESG penalties to model
ESG_environmental_penalty_constraint = m.addConstr(
    portfolio_environmental_score >= S_t_E - delta_e1 - delta_e2)
ESG_social_penalty_constraint = m.addConstr(
    portfolio_social_score >= S_t_S - delta_s1 - delta_s2)
ESG_governance_penalty_constraint = m.addConstr(
    portfolio_governance_score >= S_t_G - delta_g1 - delta_g2)
m.addConstr(delta_e1 <= 10)
m.addConstr(delta_s1 <= 10)
m.addConstr(delta_g1 <= 10);

# Add minimum return constraints
for l in scenarios:
    m.addConstr(weekly_returns.loc[l, stocks].dot(x) >= minimum_return)

    
# Add budget constraint
m.addConstr(x.sum() == 1);    

# Add blacklist constraints
m.addConstr(x.loc['ROCHE'] == 0)
m.addConstr(x.loc['NOVARTIS'] == 0);

# Add target return constraint
m.addConstr(portfolio_return >= target_return);

# Specify Gurobi options
m.setParam('OutputFlag', 0)

# Run optimization
m.optimize()

# Get portfolio return
portfolio_return.getValue()
# -

# Get weights of stocks in optimal solution
weights_p4 = pd.Series([var.X for var in x], index=stocks)
weights_p4

# Get ESG ratings P3
scores[stocks].dot(weights_p3)

# Get ESG Environmental score for P4
portfolio_environmental_score.getValue()

# Get ESG Social score for P4
portfolio_social_score.getValue()

# Get ESG Governance score for P4
portfolio_governance_score.getValue()

# ### Question:
# Assume that the model and the parameters $\alpha$ and $\mu$ are given from above. Is it possible to consturct a portfolio with its Envionmental score beeing at least 70 and 65 for the Social and Governance scores? If yes, what is the expected portfolio return?

# Adjust right-hand side of ESG penalty constraints
ESG_environmental_penalty_constraint.rhs = 70
ESG_social_penalty_constraint.rhs = 65
ESG_governance_penalty_constraint.rhs = 65

# Specify Gurobi options
m.setParam('OutputFlag', 0)

# Run optimization
m.optimize()

# Get portfolio return
portfolio_return.getValue()

# Get objective function value
m.ObjVal

# Get ESG Environmental score
portfolio_environmental_score.getValue()

# Get ESG Social score
portfolio_social_score.getValue()

# Get ESG Governance score
portfolio_governance_score.getValue()

# ### Note:
# Yes, such a portfolio is possible to construct without having to take on a penalty which is indidicated by the objective value beeing equal to zero. Further such a portfolio would have a expected return of 0.18%.

# ### Question:
# Given an investor whose main concern only consists of the highest possible Governance score? How would such a portfolio look like and what are the associated ESG-scores? 

# Adjust right-hand side of ESG penalty constraints
ESG_environmental_penalty_constraint.rhs = 0
ESG_social_penalty_constraint.rhs = 0
ESG_governance_penalty_constraint.rhs = 100

# Specify Gurobi options
m.setParam('OutputFlag', 0)

# Run optimization
m.optimize()

# Get portfolio return
portfolio_return.getValue()

# Get objective function value
m.ObjVal

# Get ESG Environmental score
portfolio_environmental_score.getValue()

# Get ESG Social score
portfolio_social_score.getValue()

# Get ESG Governance score
portfolio_governance_score.getValue()

# ### Note:
# A portfolio ESG Governance score of 73.88 can be achieved. The portfolio ESG Social score takes a value of 66.65 and the portfolio ESG Enviornmental score a value of 90.11.

# ## Perfomance evalation
# In order to evaluate our portfolios (P1-P4) and the SMI we use a second datasample containing of the weekly closing prices from the year 2018. Imagine that 1'000'000 CHF are invesed in each portflio and that it is possible to buy a fracitonal number of units.

# Import prices
weekly_prices_evaluation = pd.read_csv('prices18.csv', index_col='date')
weekly_prices_evaluation.head()

# get budget
budget=1000000

# get budget per stock
budget_per_stock_p1 = budget * weights_p1
budget_per_stock_p2 = budget * weights_p2
budget_per_stock_p3 = budget * weights_p3
budget_per_stock_p4 = budget * weights_p4

unit_per_stock_p1 = (budget_per_stock_p1 / weekly_prices_evaluation.loc['2018-01-05', stocks])
unit_per_stock_p2 = (budget_per_stock_p2 / weekly_prices_evaluation.loc['2018-01-05', stocks])
unit_per_stock_p3 = (budget_per_stock_p3 / weekly_prices_evaluation.loc['2018-01-05', stocks])
unit_per_stock_p4 = (budget_per_stock_p4 / weekly_prices_evaluation.loc['2018-01-05', stocks])
unit_smi = (budget / weekly_prices_evaluation.loc['2018-01-05', 'SMI'])

# Convert index from string to date
weekly_prices_evaluation.index = pd.to_datetime(weekly_prices_evaluation.index)

# Compute value development of portfolios
value_p1_evaluation = weekly_prices_evaluation[stocks].dot(unit_per_stock_p1)
value_p2_evaluation = weekly_prices_evaluation[stocks].dot(unit_per_stock_p2)
value_p3_evaluation = weekly_prices_evaluation[stocks].dot(unit_per_stock_p3)
value_p4_evaluation = weekly_prices_evaluation[stocks].dot(unit_per_stock_p4)
value_smi_evaluation = weekly_prices_evaluation['SMI'] * unit_smi

# %matplotlib inline
# Line plot of value development of portfolios
plt.figure(figsize=(17, 5), dpi=120)
plt.plot(value_p1_evaluation, label='P1')
plt.plot(value_p2_evaluation, label='P2')
plt.plot(value_p3_evaluation, label='P3')
plt.plot(value_p4_evaluation, label='P4')
plt.plot(value_smi_evaluation, label='SMI')
plt.axhline(1000000, color='b', linestyle='--')
plt.xlabel('Date')
plt.ylabel('Value CHF')
plt.legend();

# Get values on 28th December, 2018
print('Value P1:', value_p1_evaluation.loc['2018-12-28'].round(2))
print('Value P2:', value_p2_evaluation.loc['2018-12-28'].round(2))
print('Value P3:', value_p3_evaluation.loc['2018-12-28'].round(2))
print('Value P4:', value_p4_evaluation.loc['2018-12-28'].round(2))
print('Value SMI:', value_smi_evaluation.loc['2018-12-28'].round(2))

# ## Question:
# Compute  portfolio return and the minimum portfolio return of P1-P4 and the SMI for the years 2016-2017 as well as for the year 2018 and compare them to one another?

# ### 2016-2017

# Get unit per stock (assumption: investment of 1,000,000 CHF at beginning of year 2016)
unit_per_stock_p1 = (budget_per_stock_p1
                     / weekly_prices.loc['2016-01-01', stocks])
unit_per_stock_p2 = (budget_per_stock_p2
                     / weekly_prices.loc['2016-01-01', stocks])
unit_per_stock_p3 = (budget_per_stock_p3
                     / weekly_prices.loc['2016-01-01', stocks])
unit_per_stock_p4 = (budget_per_stock_p4
                     / weekly_prices.loc['2016-01-01', stocks])
unit_smi = (budget
            / weekly_prices.loc['2016-01-01', 'SMI'])

# Compute value development of portfolios
value_p1 = weekly_prices[stocks].dot(unit_per_stock_p1)
value_p2 = weekly_prices[stocks].dot(unit_per_stock_p2)
value_p3 = weekly_prices[stocks].dot(unit_per_stock_p3)
value_p4 = weekly_prices[stocks].dot(unit_per_stock_p4)
value_smi = weekly_prices['SMI'] * unit_smi

# Compute weekly returns for P1-P4
return_p1 = value_p1 / value_p1.shift(1) - 1
return_p2 = value_p2 / value_p2.shift(1) - 1
return_p3 = value_p3 / value_p3.shift(1) - 1
return_p4 = value_p4 / value_p4.shift(1) - 1
return_smi = value_smi / value_smi.shift(1) - 1

# Drop NA-return in first period
return_p1 = return_p1.dropna()
return_p2 = return_p2.dropna()
return_p3 = return_p3.dropna()
return_p4 = return_p4.dropna()
return_smi = return_smi.dropna()

# Get average portfolio return (in-sample)
print('In-sample')
print('Average return P1:', return_p1.mean())
print('Average return P2:', return_p2.mean())
print('Average return P3:', return_p3.mean())
print('Average return P4:', return_p4.mean())
print('Average return SMI:', return_smi.mean())

# Get minimum portfolio return (in-sample)
print('In-sample')
print('Minimum return P1:', return_p1.min())
print('Minimum return P2:', return_p2.min())
print('Minimum return P3:', return_p3.min())
print('Minimum return P4:', return_p4.min())
print('Minimum return SMI:', return_smi.min())

# ### 2018

# Compute weekly returns for P1-P4
return_p1_evaluation = value_p1_evaluation / value_p1_evaluation.shift(1) - 1
return_p2_evaluation = value_p2_evaluation / value_p2_evaluation.shift(1) - 1
return_p3_evaluation = value_p3_evaluation / value_p3_evaluation.shift(1) - 1
return_p4_evaluation = value_p4_evaluation / value_p4_evaluation.shift(1) - 1
return_smi_evaluation = value_smi_evaluation / value_smi_evaluation.shift(1) - 1

# Drop NA-return in first period
return_p1_evaluation = return_p1_evaluation.dropna()
return_p2_evaluation = return_p2_evaluation.dropna()
return_p3_evaluation = return_p3_evaluation.dropna()
return_p4_evaluation = return_p4_evaluation.dropna()
return_smi_evaluation = return_smi_evaluation.dropna()

# Get average portfolio return (out-of-sample)
print('Out-of-sample')
print('Average return P1:', return_p1_evaluation.mean())
print('Average return P2:', return_p2_evaluation.mean())
print('Average return P3:', return_p3_evaluation.mean())
print('Average return P4:', return_p4_evaluation.mean())
print('Average return SMI:', return_smi_evaluation.mean())

# Get minimum portfolio return (out-of-sample)
print('Out-of-sample')
print('Minimum return P1:', return_p1_evaluation.min())
print('Minimum return P2:', return_p2_evaluation.min())
print('Minimum return P3:', return_p3_evaluation.min())
print('Minimum return P4:', return_p4_evaluation.min())
print('Minimum return SMI:', return_smi_evaluation.min())



StarcoderdataJupyterScriptsDedupFiltered