<jupyter_start><jupyter_text>Markowitz<jupyter_code>%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt<jupyter_output><empty_output><jupyter_text>Asume que tenemos 4 acciones (assets), cada una con una serie de retorno (return) de longitud de 1000. Podemos utilizar numpy.random.randn como muestra de retorno de una distribucion normal.<jupyter_code>## NUMBER OF ASSETS
n_assets = 4

## NUMBER OF OBSERVATIONS
n_obs = 1000

return_vec = np.random.randn(n_assets, n_obs)
fig = plt.figure()
plt.plot(return_vec.T, alpha=.4);
plt.xlabel('time')
plt.ylabel('returns');<jupyter_output><empty_output><jupyter_text>Estas series de retorno pueden ser usadas para crear un amplio rango de portafolios, cual todos tienen diferentes retornos y riesgos (desviacion estandar). Podemos producir un amplio rango de vectores de pesos aleatorios y graficar esos portafolios.Queremos  que todo nuestro capital este invertido, este vector debe sumar 1.<jupyter_code>def rand_weights(n):
    ''' Produces n random weights that sum to 1 '''
    k = np.random.rand(n)
    return k / sum(k)


print(rand_weights(n_assets))
print(rand_weights(n_assets))<jupyter_output>[ 0.29889484  0.2541214   0.18900052  0.25798324]
[ 0.21954229  0.39444212  0.11827545  0.26774014]<jupyter_text>Tamaño de las acciones<jupyter_code>return_vec.shape<jupyter_output><empty_output><jupyter_text>Media de cada accion<jupyter_code>np.mean(return_vec, axis=1)
p = np.asmatrix(np.mean(return_vec, axis=1))
p<jupyter_output><empty_output><jupyter_text>Generamos un portafolio aleatorioGeneramos un vector de peso w que refleja el porcentaje de inversion en cada accion.<jupyter_code>rand_weights(return_vec.shape[0])
w = np.asmatrix(rand_weights(return_vec.shape[0]))
w<jupyter_output><empty_output><jupyter_text>Calculamos la covarianza<jupyter_code>np.cov(return_vec)
np.cov(return_vec).shape
C = np.asmatrix(np.cov(return_vec))
C<jupyter_output><empty_output><jupyter_text>El calculo del retorno con:$$R = p^Tw$$- $R$ es el retorno experado- $P^T$ es la transpuesta del vector para la media de retornos- $w$ es el vector de pesos del portafolio- $p$ es un vector columan $Nx1$ - $p^T$ es un vector fila $1xN$- $w$ es un vector columna $Nx1$ Despues calculamos el riesgo (volatilidad) con la desviacion estandar con:$$\sigma = \sqrt{w^TCw} $$- C es la matriz de covarianza de los retornos que es una matriz $NxN$**Nota:** Si calculamos la desviacion estandar con los pesos apropiados usando `std(array(ret_vec).T*w )` tendremos una pequeña diferencia en la marca `bullet`.> Esto se debe al calculo de la desviacion estandar que no toma encuenta la covarianza. >En la matriza de covarianza los valores de la diagonal representan la simple varianza de cada accion (asset) mientras que fuera de la diagonal es la varianza entre las acciones (assets)Usando la ordinaria `std()` efectivamente seria solamente la diagonal y perderiamos el resto de informacion. Pequeña per significante diferencia.<jupyter_code>np.random.seed(123)
p = np.asmatrix(np.mean(return_vec, axis=1))
w = np.asmatrix(rand_weights(return_vec.shape[0]))
C = np.asmatrix(np.cov(return_vec))
mu = w * p.T
sigma = np.sqrt(w * C * w.T)

mu,sigma
# De otra manera 
np.random.seed(123)
p = np.mean(return_vec, axis=1)
w = rand_weights(return_vec.shape[0])
C = np.cov(return_vec)
mu = p.T@w
sigma = np.sqrt(w.T @ C @ w)

mu,sigma
p.shape,w.shape,C.shape
def random_portfolio(returns):
    ''' 
    Returns the mean and standard deviation of returns for a random portfolio
    '''

    p = np.asmatrix(np.mean(returns, axis=1))
    w = np.asmatrix(rand_weights(returns.shape[0]))
    C = np.asmatrix(np.cov(returns))
    
    mu = w * p.T
    sigma = np.sqrt(w * C * w.T)
    
    # This recursion reduces outliers to keep plots pretty
    if sigma > 2:
        return random_portfolio(returns)
    return mu, sigma<jupyter_output><empty_output><jupyter_text>Vamos a generar la media de retornos y la volatilidad para 500 portafolios aleatorios.<jupyter_code>n_portfolios = 3
means, stds = np.column_stack([
    random_portfolio(return_vec) 
    for _ in range(n_portfolios)
])
means
stds
n_portfolios = 500
means, stds = np.column_stack([
    random_portfolio(return_vec) 
    for _ in range(n_portfolios)
])
fig = plt.figure()
plt.plot(stds, means, 'o', markersize=5)
plt.xlabel('std')
plt.ylabel('mean')
plt.title('Mean and standard deviation of returns of randomly generated portfolios');<jupyter_output><empty_output><jupyter_text>Una vez que tenemos la representacion de nuestros portafolios como los puntos en azules, podemos calcular la frontera de eficiencia Markowitz-style. Esto es hecho por minimizar $$w^T C w $$para $w$ en los regresos experados de los portafolios $R^T$ manteniendo la suma de todos los pesos igual a 1:$$\sum_{i}{w_i} = 1 $$Parametricamente corremos atravez de $R^Tw=μ$ y encontramos la minima varianza para diferentes $μ‘s$. Esto puede hacerse con `scipy.optimise.minimize` pero nececitamos definir un problema bastante complejo con cotas, y restricciones y multiplicadores de Lagrange.Convenientemente el paquete `cvxopt`,y un solver `convex` hace todo por nosotros.El vector `mus` produce series de valores de retornos experados $\mu$ en una nolineal(non-linear) y más apropiada maneraVermos despues que no se necesita calcular demasiados de esos y ellos se ajustan perfectamente a una parabola, cual pueden seguramente ser extramolados para valores mayores.<jupyter_code>n = len(return_vec)
returns = np.asmatrix(return_vec)
N = 100
mus = [10**(5.0 * t/N - 1.0) for t in range(N)]
plt.plot(mus)
import cvxopt as opt
from cvxopt import blas, solvers
import pandas as pd
# Turn off progress printing 
solvers.options['show_progress'] = False
np.cov(returns)
np.mean(returns, axis=1)
# Convert to cvxopt matrices
S = opt.matrix(np.cov(returns))
pbar = opt.matrix(np.mean(returns, axis=1))

S,pbar
# Create constraint matrices
G = -opt.matrix(np.eye(n))   # negative n x n identity matrix
h = opt.matrix(0.0, (n ,1))
A = opt.matrix(1.0, (1, n))
b = opt.matrix(1.0)
    
G,h,A,b<jupyter_output><empty_output><jupyter_text>mu rendimiento`p = np.mean(return_vec, axis=1)w = rand_weights(return_vec.shape[0])C = np.cov(return_vec)mu = p.T@wsigma = np.sqrt(w.T @ C @ w)` Convert to cvxopt matrices`S = opt.matrix(np.cov(returns))pbar = opt.matrix(np.mean(returns, axis=1))` p.T@w Cpbar == pS == C<jupyter_code># Calculate efficient frontier weights using quadratic programming
portfolios = [solvers.qp(mu*S, -pbar, G, h, A, b)['x'] 
                  for mu in mus]
## CALCULATE RISKS AND RETURNS FOR FRONTIER
returns = [blas.dot(pbar, x) for x in portfolios]
risks = [np.sqrt(blas.dot(x, S*x)) for x in portfolios]
## CALCULATE THE 2ND DEGREE POLYNOMIAL OF THE FRONTIER CURVE
m1 = np.polyfit(returns, risks, 2)
x1 = np.sqrt(m1[2] / m1[0])
# CALCULATE THE OPTIMAL PORTFOLIO
wt = solvers.qp(opt.matrix(x1 * S), -pbar, G, h, A, b)['x']
#np.asarray(wt), returns, risks
fig = plt.figure()
plt.plot(stds, means, 'o')
plt.ylabel('mean')
plt.xlabel('std')
plt.plot(risks, returns, 'y-o')
weights = np.asarray(wt)
print (weights)<jupyter_output>[[  4.07349921e-07]
 [  5.09503087e-01]
 [  6.77453596e-08]
 [  4.90496438e-01]]
{'fraction_non_alphanumeric': 0.11326274728336584, 'fraction_numerical': 0.026330454165505713, 'mean_word_length': 5.1622317596566525, 'pattern_counts': {'":': 0, '<': 43, '<?xml version=': 0, '>': 46, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '11', 'dir': 'jupyter-structured-clean-dedup', 'id': '128439', 'license': 'MIT', 'n_tokens_mistral': 2825, 'n_tokens_neox': 2674, 'n_words': 757, 'path': 'MarkowitzPortafolioOptimization/.ipynb_checkpoints/MarkowitzSimple-checkpoint.ipynb', 'repo_name': 'IntroCursos/MachingLearning', 'size': '7178'}