# Portfolio Optimization
> This model is an example of the classic [Markowitz portfolio selection optimization model](https://en.wikipedia.org/wiki/Markowitz_model). We want to find the fraction of the portfolio to invest among a set of stocks that balances risk and return. It is a Quadratic Programming (QP) model with vector and matrix data for returns and risk, respectively.


- toc: true 
- badges: true
- hide_binder_badge: true
- comments: true
- categories: [Investment strategy, Portfolio Optimization, Financial Modelling, Markowitz Portfolio Selection, Risk & Return, Quadratic Programming, Mathematical OPtimization, Gurobi, Python]
- image: images/when-did-you-become-an-expert-in-markowitz-portfolio-optimization-last-night.png

## Problem Statement 
The expected return and standard deviation of stock 1 are 4% and 3%, respectively. For stock 2, the expected return and standard deviation are 2% and 2%, respectively. For stock 3, the expected return and standard deviation are 1% and 1.5%, respectively. Suppose the correlation between stocks 1 and 2 is 0.2, between stocks
1 and 3 0.4, and between stocks 2 and 3 0.1.


## Model Formulation
### Sets & Indices
$i ∈ \text{stock} = \{1,2,3\}$: set of stocks
### Parameters
$\mu_i ∈ \mathbb{R}^+$:  expected return of stock $i$

$\sigma_i ∈ \mathbb{R}^+$:  standard deviation of stock $i$

$b ∈ \mathbb{R}^+$:  the minimum required expected return for the portfolio.

### Decision Variables
- $x_i$: n-element vector where each element represents the fraction of the porfolio to invest in each stock

### Objective Function
\begin{equation}
\text{Minimize} \sum_{i=1}^{n} \sum_{j=1}^{n} \mu_i \mu_j \sigma_{ij}  \tag{0}
\end{equation}

### Constraints
- **minimum expected return** 

\begin{equation}
\sum_{i=1}^{n} \mu_i x_i \geq b
\tag{1}
\end{equation}

- **total investment** 

\begin{equation}
\sum_{i=1}^{n} x_i = 1
\tag{2}
\end{equation}

- **not short-selling** 

\begin{equation}
x_i \geq 0 \quad \forall i \in \text{stock}
\tag{3}
\end{equation}



## 1) Minimum return
Find the portfolio with the lowest risk, among those with an expected return of at least 2%


```python
%pip install gurobipy
```

    Collecting gurobipy
      Downloading gurobipy-9.5.1-cp37-cp37m-manylinux2014_x86_64.whl (11.5 MB)
    [K     |████████████████████████████████| 11.5 MB 4.4 MB/s 
    [?25hInstalling collected packages: gurobipy
    Successfully installed gurobipy-9.5.1



```python
import gurobipy as gu
```


```python
# Let us start with defining the vector and matrices that contain the data.

n = 3 #number of stocks
ExpectedReturn = [0.04, 0.02, 0.01] #vector of expected returns with length n
Covariance = [[0.03, 0.00012, 0.00018], [0.00012, 0.02, 0.00003], [0.00018, 0.00003, 0.015]]
#covariance matrix with dimension nxn
MinReturn = 0.02
```


```python
# First, we create the model by function gu.Model(), where “gu” specifies that function Model() is from Gurobi
m = gu.Model("Portfolio Optimization 1")
```

    Restricted license - for non-production use only - expires 2023-10-25



```python
# Let us add the decision variables. Note that "lb" indicates the lower bound of the decision variable
StockWeight = m.addVars(n, name = 'StockWeight', lb = 0.0, vtype=gu.GRB.CONTINUOUS)
```


```python
# Let us add the objective function. Note that “for i in range(n)” indicates the loop from 0 to n-1.
m.setObjective(sum(sum(StockWeight[i]*StockWeight[j]*Covariance[i][j] for i in range(n)) for j in range(n)), gu.GRB.MINIMIZE)
```


```python
# Let us add the constaints.
m.addConstr(sum(StockWeight[i] for i in range(n)) == 1, name = 'SumofWeights')
m.addConstr(sum(StockWeight[i]*ExpectedReturn[i] for i in range(n)) >= MinReturn, name = 'ConstraintonMinReturn')
```




    <gurobi.Constr *Awaiting Model Update*>




```python
m.optimize()
```

    Gurobi Optimizer version 9.5.1 build v9.5.1rc2 (linux64)
    Thread count: 1 physical cores, 2 logical processors, using up to 2 threads
    Optimize a model with 2 rows, 3 columns and 6 nonzeros
    Model fingerprint: 0x390fa4c2
    Model has 6 quadratic objective terms
    Coefficient statistics:
      Matrix range     [1e-02, 1e+00]
      Objective range  [0e+00, 0e+00]
      QObjective range [1e-04, 6e-02]
      Bounds range     [0e+00, 0e+00]
      RHS range        [2e-02, 1e+00]
    Presolve time: 0.02s
    Presolved: 2 rows, 3 columns, 6 nonzeros
    Presolved model has 6 quadratic objective terms
    Ordering time: 0.00s
    
    Barrier statistics:
     Free vars  : 2
     AA' NZ     : 6.000e+00
     Factor NZ  : 1.000e+01
     Factor Ops : 3.000e+01 (less than 1 second per iteration)
     Threads    : 1
    
                      Objective                Residual
    Iter       Primal          Dual         Primal    Dual     Compl     Time
       0   3.20152377e+05 -3.20152377e+05  2.50e+03 2.47e-07  1.00e+06     0s
       1   1.15177644e+04 -1.17907196e+04  1.97e+02 1.95e-08  8.42e+04     0s
       2   4.34195814e-02 -3.58271785e+02  6.80e-01 6.72e-11  3.82e+02     0s
       3   6.98591118e-03 -1.98336951e+02  6.80e-07 6.94e-17  4.96e+01     0s
       4   6.98585719e-03 -1.92714187e-01  4.66e-12 0.00e+00  4.99e-02     0s
       5   6.95633171e-03  5.66424294e-03  1.48e-14 0.00e+00  3.23e-04     0s
       6   6.75813420e-03  6.68576945e-03  6.94e-17 6.94e-18  1.81e-05     0s
       7   6.73329099e-03  6.72497295e-03  3.41e-15 2.74e-18  2.08e-06     0s
       8   6.72958274e-03  6.72842709e-03  1.83e-14 1.25e-18  2.89e-07     0s
       9   6.72896822e-03  6.72882136e-03  5.30e-14 6.94e-18  3.67e-08     0s
      10   6.72886123e-03  6.72885109e-03  6.70e-14 2.74e-19  2.54e-09     0s
      11   6.72885183e-03  6.72885180e-03  1.17e-13 6.94e-18  6.72e-12     0s
    
    Barrier solved model in 11 iterations and 0.05 seconds (0.00 work units)
    Optimal objective 6.72885183e-03
    



```python
# Let us have a look at the solution.
var = m.getVars()
for i in range(n):
  print(var[i].varName, var[i].x)
```

    StockWeight[0] 0.221851545261853
    StockWeight[1] 0.3344459277157383
    StockWeight[2] 0.4437025270221763



```python
print('Objective function', m.objVal)
```

    Objective function 0.00672885182871444


## 2) What if we cannot invest more than half of our money in stocks 1 and 2?
To answer part 2, we need to add a new constraint involving stocks’ 1 and 2 decision variables.

### Parameters
$a ∈ \mathbb{R}^+$:  the maximum required investment for the stock 1 and 2.

### Constraints
- **Maximum investment** 

\begin{equation}
x_1 + x_2 = a
\tag{1}
\end{equation}



```python
# Let us start with defining the vector and matrices that contain the data.
MaxInvest = 0.5
```


```python
# Let us add the constaints.
m.addConstr(StockWeight[0] + StockWeight[1] == MaxInvest, name = 'MaxInvest')
```




    <gurobi.Constr *Awaiting Model Update*>




```python
m.optimize()
```

    Gurobi Optimizer version 9.5.1 build v9.5.1rc2 (linux64)
    Thread count: 1 physical cores, 2 logical processors, using up to 2 threads
    Optimize a model with 3 rows, 3 columns and 8 nonzeros
    Model fingerprint: 0xb98dfd3d
    Model has 6 quadratic objective terms
    Coefficient statistics:
      Matrix range     [1e-02, 1e+00]
      Objective range  [0e+00, 0e+00]
      QObjective range [1e-04, 6e-02]
      Bounds range     [0e+00, 0e+00]
      RHS range        [2e-02, 1e+00]
    Presolve time: 0.02s
    Presolved: 3 rows, 3 columns, 8 nonzeros
    Presolved model has 6 quadratic objective terms
    Ordering time: 0.00s
    
    Barrier statistics:
     Free vars  : 2
     AA' NZ     : 1.000e+01
     Factor NZ  : 1.500e+01
     Factor Ops : 5.500e+01 (less than 1 second per iteration)
     Threads    : 1
    
                      Objective                Residual
    Iter       Primal          Dual         Primal    Dual     Compl     Time
       0   3.20072058e+05 -3.20072058e+05  2.50e+03 6.01e-02  1.00e+06     0s
       1   4.83954583e+03 -5.08987963e+03  2.54e+02 6.11e-03  1.05e+05     0s
       2   3.32693448e-02 -3.05887397e+02  5.42e-01 1.30e-05  3.01e+02     0s
       3   6.97658803e-03 -1.15906961e+02  5.42e-07 1.30e-11  2.90e+01     0s
       4   6.97658263e-03 -1.09052744e-01  5.43e-13 1.35e-14  2.90e-02     0s
       5   6.97397122e-03  6.75207038e-03  1.11e-15 1.04e-16  5.55e-05     0s
       6   6.94299555e-03  6.94200376e-03  1.11e-16 6.94e-18  2.48e-07     0s
       7   6.94250050e-03  6.94249950e-03  8.33e-16 6.94e-18  2.48e-10     0s
    
    Barrier solved model in 7 iterations and 0.05 seconds (0.00 work units)
    Optimal objective 6.94250050e-03
    



```python
# Let us have a look at the solution.
var = m.getVars()
for i in range(n):
  print(var[i].varName, var[i].x)
```

    StockWeight[0] 0.25000009613332785
    StockWeight[1] 0.24999990386667328
    StockWeight[2] 0.5000000000000022



```python
print('Objective function', m.objVal)
```

    Objective function 0.006942500495087108


c11c1b649614f4f433b0e8de6d37ba30b89c9613
ipynb
Jupyter Notebook
_notebooks/2022_05_02.ipynb
TheMuahDibPortfolio/Blog_Finance
51259216f0be90345e5629757cecb2752024bdfb
Apache-2.0
_notebooks/2022_05_02.ipynb
TheMuahDibPortfolio/Blog_Finance
51259216f0be90345e5629757cecb2752024bdfb
Apache-2.0
_notebooks/2022_05_02.ipynb
TheMuahDibPortfolio/Blog_Finance
51259216f0be90345e5629757cecb2752024bdfb
Apache-2.0
Qwen/Qwen-72B
1. YES
2. YES