From 7bda300ec4b9402480983ddcb0580223c1cf93f7 Mon Sep 17 00:00:00 2001 From: Cursor Agent Date: Sat, 20 Jun 2026 01:35:09 +0000 Subject: [PATCH] Fix cardinality notebook execution Co-authored-by: Antonije Mirkovic --- Cardinality_Constrained.ipynb | 1060 ++++++++++++++++++--------------- README.md | 7 +- requirements.txt | 3 + 3 files changed, 579 insertions(+), 491 deletions(-) create mode 100644 requirements.txt diff --git a/Cardinality_Constrained.ipynb b/Cardinality_Constrained.ipynb index d108cd2..be8786c 100644 --- a/Cardinality_Constrained.ipynb +++ b/Cardinality_Constrained.ipynb @@ -1,490 +1,574 @@ { - "cells": [ - { - "cell_type": "markdown", - "id": "c0845ac3", - "metadata": {}, - "source": [ - "# Cardinality Constrained Portfolio Optimization\n", - "\n", - "This notebook explores the **Cardinality-Constrained Mean-Variance (CCMV)** portfolio optimization problem, based on the methodology presented in *\"Cardinality-Constrained Portfolio Optimization with Clustering\"* (Ebrahimi et al., 2025).\n", - "\n", - "The classical **Markowitz mean-variance model** is extended by adding a **cardinality constraint**, limiting the number of assets included in the portfolio. This transforms the problem into a **Mixed-Integer Quadratic Program (MIQP)**, making it significantly more challenging to solve due to its NP-hard nature.\n", - "\n", - "In this notebook, we will:\n", - "\n", - "- Formulate the classical and cardinality-constrained optimization problems\n", - "- Discuss the implications of binary variables and computational complexity\n", - "- Implement and solve small CCMV problems using Python\n", - "- Prepare for future integration of clustering methods as described in the paper\n", - "\n", - "This work serves as both a theoretical deep dive and a practical foundation for reproducing the paper’s results." - ] + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Cardinality Constrained Portfolio Optimization\n", + "\n", + "This notebook explores the **Cardinality-Constrained Mean-Variance (CCMV)** portfolio optimization problem, based on the methodology presented in *\"Cardinality-Constrained Portfolio Optimization with Clustering\"* (Ebrahimi et al., 2025).\n", + "\n", + "The classical **Markowitz mean-variance model** is extended by adding a **cardinality constraint**, limiting the number of assets included in the portfolio. This transforms the problem into a **Mixed-Integer Quadratic Program (MIQP)**, making it significantly more challenging to solve due to its NP-hard nature.\n", + "\n", + "In this notebook, we will:\n", + "\n", + "- Formulate the classical and cardinality-constrained optimization problems\n", + "- Discuss the implications of binary variables and computational complexity\n", + "- Implement and solve small CCMV problems using Python\n", + "- Prepare for future integration of clustering methods as described in the paper\n", + "\n", + "This work serves as both a theoretical deep dive and a practical foundation for reproducing the paper’s results." + ], + "id": "c0845ac3" + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Portfolio Optimization Problem Formulation\n", + "\n", + "The classical **mean-variance portfolio optimization** problem, introduced by Markowitz (1952), aims to minimize portfolio risk for a given level of expected return (or equivalently, balance both using a risk aversion parameter). It is formulated as:\n", + "\n", + "$$\n", + "\\min_{w \\in \\mathbb{R}^N} \\quad w^\\top \\Sigma w - \\gamma \\mu^\\top w \\quad \\text{subject to} \\quad \\mathbf{1}^\\top w = 1,\\quad 0 \\leq w_i \\leq 1 \\quad \\forall i \\in \\{1, \\ldots, N\\}\n", + "$$\n", + "\n", + "where:\n", + "\n", + "- $ w = (w_1, \\ldots, w_N)^\\top $ is the vector of portfolio weights \n", + "- $ \\Sigma \\in \\mathbb{R}^{N \\times N} $ is the **covariance matrix** of asset returns \n", + "- $ \\mu \\in \\mathbb{R}^N $ is the **expected return vector** \n", + "- $ \\gamma > 0 $ is the **risk aversion coefficient** \n", + "- $ \\mathbf{1} \\in \\mathbb{R}^N $ is a vector of ones \n", + "- The constraint $ \\mathbf{1}^\\top w = 1 $ ensures the portfolio is fully invested \n", + "- The box constraint $ 0 \\leq w_i \\leq 1 $ limits individual asset weights\n", + "\n", + "This is a **convex quadratic program** and can be efficiently solved using standard optimization methods.\n", + "\n" + ], + "id": "71b7d252" + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "This classical model is easy to solve due to its convexity, but it does not capture important practical constraints that investors often face. In real-world applications, it is common to:\n", + "\n", + "- Reduce transaction costs \n", + "- Control portfolio complexity \n", + "- Avoid excessive concentration in a few assets\n", + "\n", + "To reflect these needs, we introduce two key extensions:\n", + "\n", + "#### 1. Cardinality Constraint\n", + "\n", + "Investors may want to limit the total number of assets in the portfolio. Let $\\Delta \\in \\mathbb{N}_+$ denote the maximum allowed number of assets. This leads to the **cardinality constraint**:\n", + "\n", + "$$\n", + "\\|w\\|_0 \\leq \\Delta\n", + "$$\n", + "\n", + "Here, $\\|w\\|_0$ is the **zero-norm**, which counts the number of nonzero entries in $w$. This constraint is **non-convex** and introduces **combinatorial complexity**, requiring the use of binary decision variables.\n", + "\n", + "#### 2. Box Constraint on Asset Weights\n", + "\n", + "To prevent overexposure to individual assets, we constrain each weight to lie within a maximum threshold. Let $\\bar{w} \\in (0, 1]$ be the maximum allowable weight per asset. The constraint is written as:\n", + "\n", + "$$\n", + "w \\in [0, \\bar{w}]^N\n", + "$$\n", + "\n", + "This **box constraint** ensures no single asset can dominate the portfolio, and also helps with regulatory compliance or risk diversification.\n", + "\n", + "Together, these constraints transform the problem into a **Mixed-Integer Quadratic Program (MIQP)**, which is computationally much harder to solve than the unconstrained case.\n", + "\n", + "We can then formulate the **Cardinality-Constrained Mean-Variance (CCMV)** optimization problem as follows:\n", + "\n", + "$$\n", + "w^*(\\Delta, \\bar{w}) = \\arg\\min_{w \\in \\mathbb{R}^N} \\quad w^\\top \\Sigma w - \\gamma \\mu^\\top w \\quad \\text{subject to} \\quad \\mathbf{1}^\\top w = 1, \\quad 0 \\leq w_i \\leq \\bar{w} \\quad \\forall i, \\quad \\|w\\|_0 \\leq \\Delta\n", + "$$\n", + "\n", + "Furthermore, it can be shown that this problem is NP-hard due to the combinatorial nature of the cardinality constraint. The introduction of binary variables to enforce this constraint complicates the optimization significantly.\n" + ], + "id": "a7602ae4" + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Solving the CCMV Problem\n", + "\n", + "The **Cardinality-Constrained Mean-Variance (CCMV)** problem is computationally challenging due to its non-convex nature and the presence of binary decisions. To address this, we adopt the **clustering-based approach** proposed in the paper.\n", + "\n", + "Specifically, we apply **hierarchical clustering** to group assets based on the similarity of their return profiles. Each cluster represents a group of highly correlated assets, and we solve the CCMV problem using the **cluster centroids** as representative assets.\n", + "\n", + "This technique offers two key advantages:\n", + "- It **reduces the dimensionality** of the problem, making it more tractable.\n", + "- It **preserves key structure** in the data by exploiting the correlation between assets within each cluster.\n", + "\n", + "By solving a reduced CCMV problem on the cluster level, we approximate the original high-dimensional problem while significantly improving computational efficiency.\n" + ], + "id": "de9a9d14" + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Implementation of Hierarchical Clustering \n", + "\n", + "We now implement the hierarchical clustering approach, following the pseudo-code provided in the paper. This algorithm clusters assets based on the similarity of their return profiles by analyzing the correlation structure of the return matrix.\n", + "\n", + "The purpose of this step is to reduce the dimensionality of the portfolio optimization problem by grouping similar assets together. This allows us to simplify the later stages of solving the CCMV problem, while still preserving the core structure of the asset universe.\n", + "\n", + "The clustering is performed using a bottom-up (agglomerative) approach, with distances computed between cluster centroids in the distance space. The final partition is determined by cutting the dendrogram at the optimal level, based on the largest gap in merge distances." + ], + "id": "cbf7e6c1" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "import itertools\n", + "\n", + "import numpy as np\n", + "from scipy.optimize import brentq, minimize\n", + "import matplotlib.pyplot as plt" + ], + "execution_count": 2, + "outputs": [], + "id": "95332e05" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "def Hierarchical_Clustering(R: np.ndarray) -> list[list[int]]:\n", + "\t\"\"\"\n", + "\tPerform agglomerative clustering on the given return matrix R.\n", + "\t\n", + "\tParameters:\n", + "\tR (np.ndarray): A return matrix with observations in rows and assets in columns.\n", + "\t\n", + "\tReturns:\n", + "\tlist[list[int]]: Asset-index clusters.\n", + "\t\"\"\"\n", + "\treturns = np.asarray(R)\n", + "\tif returns.ndim != 2:\n", + "\t\traise ValueError(\"R must be a two-dimensional return matrix\")\n", + "\n", + "\tn_assets = returns.shape[1]\n", + "\tif n_assets == 0:\n", + "\t\treturn []\n", + "\tif n_assets == 1:\n", + "\t\treturn [[0]]\n", + "\n", + "\t# 1. Compute the correlation matrix tau. Constant columns can produce NaNs.\n", + "\ttau = np.nan_to_num(np.corrcoef(returns, rowvar=False), nan=0.0)\n", + "\tnp.fill_diagonal(tau, 1.0)\n", + "\n", + "\t# 2. Construct the distance matrix D\n", + "\tD = 1 - tau\n", + "\n", + "\t# 3. Initialize clusters C\n", + "\tC = [[i] for i in range(n_assets)]\n", + "\n", + "\tdef d(Ci, Cj):\n", + "\t\t\"\"\"\n", + "\t\tCompute Ward-style distance between two clusters Ci and Cj.\n", + "\t\t\"\"\"\n", + "\t\tcent_i = np.sum(D[Ci, :], axis=0) / len(Ci)\n", + "\t\tcent_j = np.sum(D[Cj, :], axis=0) / len(Cj)\n", + "\t\tdist = (len(Ci) * len(Cj)) / (len(Ci) + len(Cj)) * np.linalg.norm(cent_i - cent_j) ** 2\n", + "\t\treturn dist\n", + "\n", + "\t# 4. Repeat until only one cluster remains\n", + "\tmerge_distances = []\n", + "\tmerge_log = []\n", + "\twhile len(C) > 1:\n", + "\t\tmin_distance = float('inf')\n", + "\t\tbest_pair = (0, 0)\n", + "\t\tfor i in range(len(C)):\n", + "\t\t\tfor j in range(i + 1, len(C)):\n", + "\t\t\t\tdistance = d(C[i], C[j])\n", + "\t\t\t\tif distance < min_distance:\n", + "\t\t\t\t\tmin_distance = distance\n", + "\t\t\t\t\tbest_pair = (i, j)\n", + "\t\t\n", + "\t\t# Merge the best pair of clusters\n", + "\t\ti, j = best_pair\n", + "\t\tnew_cluster = C[i] + C[j]\n", + "\t\tmerge_log.append((C[i], C[j], new_cluster, min_distance))\n", + "\t\tdel C[max(i, j)] # Remove the larger index first to avoid index issues\n", + "\t\tdel C[min(i, j)]\n", + "\t\tC.append(new_cluster)\n", + "\t\tmerge_distances.append(min_distance)\n", + "\n", + "\t# 5. Determine the cut point by stopping before the largest jump in merge distance.\n", + "\tif len(merge_distances) == 1:\n", + "\t\tthreshold = merge_distances[0]\n", + "\telse:\n", + "\t\tlargest_gap_idx = int(np.argmax(np.diff(merge_distances)))\n", + "\t\tthreshold = merge_distances[largest_gap_idx]\n", + "\n", + "\t# 6. Replay merges up to the threshold to recover asset-index clusters.\n", + "\tclusters = [{i} for i in range(n_assets)]\n", + "\tfor left, right, _, distance in merge_log:\n", + "\t\tif distance > threshold:\n", + "\t\t\tcontinue\n", + "\n", + "\t\tleft_set = set(left)\n", + "\t\tright_set = set(right)\n", + "\t\tto_merge = [\n", + "\t\t\tcluster\n", + "\t\t\tfor cluster in clusters\n", + "\t\t\tif cluster & left_set or cluster & right_set\n", + "\t\t]\n", + "\t\tif len(to_merge) < 2:\n", + "\t\t\tcontinue\n", + "\n", + "\t\tmerged = set().union(*to_merge)\n", + "\t\tfor cluster in to_merge:\n", + "\t\t\tclusters.remove(cluster)\n", + "\t\tclusters.append(merged)\n", + "\n", + "\treturn [sorted(cluster) for cluster in clusters]" + ], + "execution_count": 3, + "outputs": [], + "id": "88e623bf" + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Modeling the Clustered CCMV Problem\n", + "\n", + "Let $C = \\{C_1, \\dots, C_K\\}$ be the optimal clusters obtained from hierarchical clustering. We now introduce three key parameters:\n", + "\n", + "- Maximum allowable weight for any asset $\\bar{w}$\n", + "- Maximum number of assets per cluster $\\Delta = (\\Delta_1, \\dots \\Delta_k)$\n", + "- $\\alpha = (\\alpha_1, \\dots, \\alpha_k)$, where $\\alpha_i$ is the fraction of the portfolio allocated to cluster $C_i$\n", + "\n", + "We can then define the feasible set of portfolio weights accounting for both clustering and cardinality constraints:\n", + "$$\n", + "\t\\mathcal{A}_C(\\Delta, \\alpha, \\bar{w}) = \\left\\{ w \\in [0, \\bar{w}]^N : \\|w_{C_i}\\|_0 \\leq \\Delta_i, \\quad \\mathbf{1}^\\top w_{C_i} = \\alpha_i \\quad \\forall i = 1, \\ldots, K \\right\\}\n", + "\n", + "$$\n", + "\n", + "We can then express the clustered CCMV problem as:\n", + "\n", + "$$\n", + "\tw^*(C, \\Delta, \\alpha, \\bar{w}) = \\arg\\min_{w \\in \\mathcal{A}_C(\\Delta, \\alpha, \\bar{w})} \\quad w^\\top \\Sigma w - \\gamma \\mu^\\top w\n", + "$$\n", + "\n", + "Our goal is to determine the optimal cluster-level allocation $\\alpha_k$, the maximum number of assets $\\Delta_k$ per cluster $C_k$, and the optimal asset weights $w$ to balance allocation and enchance portfolio diversification.\n", + "\n", + "We now introduce an intra-cluster optimization problem to determine the optimal asset weights within each cluster. This is formulated as:\n", + "$$\n", + "\tw^*_{C_k} = \\arg\\min_{w \\in \\mathcal{A}_{C_k}(\\Delta_k, \\alpha_k, \\bar{w})} \\quad w^\\top \\Sigma_{C_k} w - \\gamma \\mu_{C_k}^\\top w\n", + "$$\n", + "\n", + "Where $\\Sigma_{C_k}$ and $\\mu_{C_k}$ are the covariance matrix and expected return vector for the assets in cluster $C_k$. We can also define the feasible set used above as:\n", + "$$\n", + "\\mathcal{A}_{C_k}(\\Delta_k, \\alpha_k, \\bar{w}) = \\left\\{ w \\in [0, \\bar{w}]^{|C_k|} : \\|w\\|_0 \\leq \\Delta_k, \\quad 1^Tw = \\alpha_k \\right\\}.\n", + "$$\n", + "\n", + "This means that we can write the original feasible set as:\n", + "$$\n", + "\\mathcal{A}_C(\\Delta, \\alpha, \\bar{w}) = \\bigcap_{k=1}^K \\left\\{ w : w_{C_k} \\in \\mathcal{A}_{C_k}(\\Delta_k, \\alpha_k, \\bar{w}) \\right\\}.\n", + "$$\n", + "\n", + "We now look at the objective function. To simplify things, we neglect correlations between different clusters, which allows us to write the objective function as a sum of intra-cluster objectives:\n", + "$$\n", + "f(w) \\approx \\sum_{k=1}^K f_{C_k}(w) = \\sum_{k=1}^K \\left( \\sum_{i,j \\in C_k} w_iw_j\\sigma_{ij} -\\gamma \\sum_{i \\in C_k} w_i \\mu_i \\right).\n", + "$$\n", + "\n", + "We can then approximate the optimal objective value as \n", + "$$\n", + "f^*_C(\\Delta, \\alpha) \\approx \\sum_{k=1}^K f^*_{C_k}(\\Delta_k, \\alpha_k)\n", + "$$\n", + "\n", + "where $f^*_C(\\Delta, \\alpha) = \\min_{w \\in \\mathcal{A}_C(\\Delta, \\alpha, \\bar{w})} f(w)$ and $f^*_{C_k}(\\Delta_k, \\alpha_k) = \\min_{w_{C_k} \\in \\mathcal{A}_{C_k}(\\Delta_k, \\alpha_k, \\bar{w})} f_{C_k}(w)$.\n", + "\n", + "We can then finally define\n", + "\n", + "$$\n", + "\tf^*_C(\\Delta) := \\min_{\\substack{\\Delta, \\alpha, w_i; \\\\ \\sum_{k=1}^K \\Delta_k = \\Delta; \\\\ \\sum_{k=1}^K \\alpha_k = 1}} f^*_C(\\Delta, \\alpha) = \\min_{\\substack{\\Delta, \\alpha, w_i; \\\\ \\sum_{k=1}^K \\Delta_k = \\Delta; \\\\ \\sum_{k=1}^K \\alpha_k = 1}} \\min_{w \\in \\mathcal{A}_C(\\Delta, \\alpha, \\bar{w})} f(w)\n", + "$$\n", + "\n", + "This final formulation represents a two-level optimization: the outer minimization selects how to allocate weight ($\\alpha$) and sparsity ($\\Delta$) across clusters, while the inner minimization solves for the optimal asset weights within each cluster under those constraints. Together, this balances diversification across clusters with precise selection within them.\n" + ], + "id": "40b512f4" + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Implementation of $\\Delta$-CCMV and $\\alpha$-CCMV\n", + "\n", + "We now implement the $\\Delta$-CCMV and $\\alpha$-CCMV problems as described by the pseuducode in the paper. The first method $\\Delta$-CCMV pre-assigns the number of selected assets per cluster, while the second method $\\alpha$-CCMV pre-assigns the fraction of the portfolio allocated to each cluster." + ], + "id": "dab701ff" + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### $\\Delta$-CCMV\n", + "\n", + "We introduce a score based method to determine the number of selected assets per cluster. For each cluster $C_k$, we compute a score based on the expected return and risk of the assets within that cluster.\n", + "We compute a continuous allocation $\\Delta \\times s_k$ for each cluster $C_k$, where $s_k$ is the relative importance of the cluster, then convert those values to integer cardinalities with a largest-remainder rule so the final allocations sum to $\\Delta$. In the algorithm, we compute the score as follows:\n", + "$$\n", + "s_k = \\left(\\frac{\\gamma \\mu_{C_k} - \\lambda}{2 \\sigma^2_{C_k}}\\right)_+\n", + "$$\n", + "\n", + "where $\\mu_{C_k}$ and $\\sigma^2_{C_k}$ are the expected return and variance of the assets in cluster $C_k$. $\\lambda$ is a normalization constant that ensures the scores sum to 1.\n", + "\n", + "By using scores, we can ensure that the clusters with stronger risk-adjusted performance are more represented in the final portfolio, while still respecting the cardinality constraint." + ], + "id": "d751606d" + }, + { + "cell_type": "code", + "metadata": {}, + "source": [ + "def allocate_cardinality(clusters, mu, Sigma, Delta, gamma):\n", + " \"\"\"\n", + " Allocate the number of selected assets (Delta_k) to each cluster.\n", + "\n", + " The continuous cluster scores sum to one. They are converted to integer\n", + " cardinalities with a largest-remainder rule while respecting cluster sizes.\n", + " \"\"\"\n", + " mu = np.asarray(mu)\n", + " Sigma = np.asarray(Sigma)\n", + " Delta = int(Delta)\n", + "\n", + " if Delta < 0:\n", + " raise ValueError(\"Delta must be non-negative\")\n", + " if len(clusters) == 0:\n", + " if Delta == 0:\n", + " return []\n", + " raise ValueError(\"Cannot allocate cardinality without clusters\")\n", + "\n", + " capacities = np.array([len(cluster) for cluster in clusters], dtype=int)\n", + " total_capacity = int(capacities.sum())\n", + " if Delta > total_capacity:\n", + " raise ValueError(\"Delta cannot exceed the number of available assets\")\n", + " if Delta == 0:\n", + " return [0 for _ in clusters]\n", + "\n", + " # Step 1: Compute cluster-level expected returns and variances\n", + " cluster_returns = []\n", + " cluster_variances = []\n", + " for cluster in clusters:\n", + " cluster_assets = np.array(cluster)\n", + " cluster_mu = mu[cluster_assets]\n", + " cluster_Sigma = Sigma[np.ix_(cluster_assets, cluster_assets)]\n", + " \n", + " expected_return = float(np.mean(cluster_mu))\n", + " variance = max(float(np.mean(np.diag(cluster_Sigma))), 1e-12)\n", + " \n", + " cluster_returns.append(expected_return)\n", + " cluster_variances.append(variance)\n", + " \n", + " # Step 2: Define the score function\n", + " def score_function(lambda_, mu_k, sigma2_k, gamma):\n", + " scores = []\n", + " for mu_i, sigma2_i in zip(mu_k, sigma2_k):\n", + " raw = (gamma * mu_i - lambda_) / (2 * sigma2_i)\n", + " scores.append(max(0, raw))\n", + " return np.array(scores, dtype=float)\n", + "\n", + " def root_objective(lambda_, mu_k, sigma2_k, gamma):\n", + " return float(score_function(lambda_, mu_k, sigma2_k, gamma).sum() - 1)\n", + "\n", + " # Step 3: Solve for lambda*. Expand the bracket from the data scale.\n", + " center = float(np.max(gamma * np.array(cluster_returns)))\n", + " span = max(1.0, abs(center), float(np.max(cluster_variances)))\n", + " lower = center - span\n", + " upper = center + span\n", + " while root_objective(lower, cluster_returns, cluster_variances, gamma) <= 0:\n", + " lower -= span\n", + " span *= 2\n", + " while root_objective(upper, cluster_returns, cluster_variances, gamma) >= 0:\n", + " upper += span\n", + " span *= 2\n", + "\n", + " lambda_star = brentq(root_objective, lower, upper, args=(cluster_returns, cluster_variances, gamma))\n", + "\n", + " # Step 4: Final scores and integer allocation\n", + " scores = score_function(lambda_star, cluster_returns, cluster_variances, gamma)\n", + " raw_allocations = Delta * scores\n", + " allocations = np.minimum(np.floor(raw_allocations).astype(int), capacities)\n", + "\n", + " remaining = Delta - int(allocations.sum())\n", + " remainders = raw_allocations - np.floor(raw_allocations)\n", + " while remaining > 0:\n", + " eligible = np.where(allocations < capacities)[0]\n", + " if len(eligible) == 0:\n", + " raise ValueError(\"Unable to allocate all cardinality within cluster capacities\")\n", + " best_idx = max(eligible, key=lambda idx: (remainders[idx], scores[idx], capacities[idx]))\n", + " allocations[best_idx] += 1\n", + " remaining -= 1\n", + "\n", + " return allocations.astype(int).tolist()\n", + "\n", + "\n", + "def solve_intracluster_miqp(mu_cluster, Sigma_cluster, delta_k, alpha_k, bar_w, gamma):\n", + " \"\"\"\n", + " Solve a cardinality-constrained intra-cluster problem by enumerating subsets.\n", + "\n", + " For the small examples targeted by this notebook, enumerating candidate asset\n", + " subsets and solving the continuous convex subproblem with SciPy avoids relying\n", + " on an external commercial MIQP solver.\n", + " \"\"\"\n", + " mu_cluster = np.asarray(mu_cluster)\n", + " Sigma_cluster = np.asarray(Sigma_cluster)\n", + " n = len(mu_cluster)\n", + " delta_k = int(delta_k)\n", + " alpha_k = float(alpha_k)\n", + "\n", + " if delta_k < 0:\n", + " raise ValueError(\"delta_k must be non-negative\")\n", + " if alpha_k < -1e-10:\n", + " raise ValueError(\"alpha_k must be non-negative\")\n", + " if alpha_k <= 1e-10:\n", + " return np.zeros(n)\n", + " if delta_k == 0:\n", + " raise ValueError(\"Positive cluster allocation requires positive delta_k\")\n", + " if alpha_k > delta_k * bar_w + 1e-10:\n", + " raise ValueError(\"Cluster allocation exceeds the selected assets' weight capacity\")\n", + "\n", + " best_value = float('inf')\n", + " best_weights = None\n", + " max_subset_size = min(delta_k, n)\n", + "\n", + " for subset_size in range(1, max_subset_size + 1):\n", + " if alpha_k > subset_size * bar_w + 1e-10:\n", + " continue\n", + "\n", + " for subset in itertools.combinations(range(n), subset_size):\n", + " subset = np.array(subset)\n", + " Sigma_sub = Sigma_cluster[np.ix_(subset, subset)]\n", + " mu_sub = mu_cluster[subset]\n", + " x0 = np.full(subset_size, alpha_k / subset_size)\n", + "\n", + " def objective(x):\n", + " return float(x @ Sigma_sub @ x - gamma * mu_sub @ x)\n", + "\n", + " def gradient(x):\n", + " return 2 * Sigma_sub @ x - gamma * mu_sub\n", + "\n", + " constraints = ({\n", + " 'type': 'eq',\n", + " 'fun': lambda x: float(np.sum(x) - alpha_k),\n", + " 'jac': lambda x: np.ones_like(x),\n", + " },)\n", + " bounds = [(0.0, bar_w) for _ in range(subset_size)]\n", + " result = minimize(\n", + " objective,\n", + " x0,\n", + " jac=gradient,\n", + " bounds=bounds,\n", + " constraints=constraints,\n", + " method='SLSQP',\n", + " options={'ftol': 1e-10, 'maxiter': 500},\n", + " )\n", + " if not result.success:\n", + " continue\n", + "\n", + " value = objective(result.x)\n", + " if value < best_value:\n", + " weights = np.zeros(n)\n", + " weights[subset] = np.where(np.abs(result.x) < 1e-10, 0.0, result.x)\n", + " best_value = value\n", + " best_weights = weights\n", + "\n", + " if best_weights is None:\n", + " raise RuntimeError(\"No feasible intra-cluster portfolio found\")\n", + "\n", + " return best_weights\n", + "\t\n", + "\n", + "def delta_ccmv(R, mu, Sigma, Delta, bar_w, gamma):\n", + "\t\"\"\"\n", + "\tSolve the CCMV problem by pre-assigning the number of assets in each cluster.\n", + "\t\n", + "\tParameters:\n", + "\tR (np.ndarray): A return matrix.\n", + "\tmu (np.ndarray): Expected returns.\n", + "\tSigma (np.ndarray): Covariance matrix of returns.\n", + "\tDelta (int): Total number of assets to include in the portfolio (cardinality constraint)\n", + "\tbar_w (float): Maximum allowed weight per asset\n", + "\tgamma (float): Risk aversion parameter.\n", + "\t\n", + "\tReturns:\n", + "\tw: Optimal portfolio weights\n", + "\talpha: cluster allocation proportions\n", + "\t\"\"\"\n", + "\tR = np.asarray(R)\n", + "\tmu = np.asarray(mu)\n", + "\tSigma = np.asarray(Sigma)\n", + "\n", + "\tif R.ndim != 2:\n", + "\t\traise ValueError(\"R must be a two-dimensional return matrix\")\n", + "\tif R.shape[1] != len(mu) or Sigma.shape != (len(mu), len(mu)):\n", + "\t\traise ValueError(\"R, mu, and Sigma dimensions are inconsistent\")\n", + "\tif Delta < 1:\n", + "\t\traise ValueError(\"Delta must be at least 1 for a fully invested portfolio\")\n", + "\tif Delta * bar_w < 1 - 1e-10:\n", + "\t\traise ValueError(\"Delta * bar_w must be at least 1 to invest the full portfolio\")\n", + "\n", + "\t# Step 1: Collect the clusters from R\n", + "\tclusters = Hierarchical_Clustering(R)\n", + "\n", + "\t# Step 2: Delta pre-assignment\n", + "\tDelta_k = allocate_cardinality(clusters, mu, Sigma, Delta, gamma)\n", + " \n", + "\t# Step 3: Allocate portfolio weight in proportion to each cluster's capacity.\n", + "\tcluster_capacities = np.array(Delta_k, dtype=float) * bar_w\n", + "\tif cluster_capacities.sum() < 1 - 1e-10:\n", + "\t\traise ValueError(\"Allocated cluster capacity is insufficient for a fully invested portfolio\")\n", + "\talpha_k = (cluster_capacities / cluster_capacities.sum()).tolist()\n", + "\n", + "\tw = np.zeros(R.shape[1])\n", + "\tfor i, cluster in enumerate(clusters):\n", + "\t\tmu_cluster = mu[cluster]\n", + "\t\tSigma_cluster = Sigma[np.ix_(cluster, cluster)]\n", + "\t\tdelta_k = Delta_k[i]\n", + "\t\talpha_k_val = alpha_k[i]\n", + "\t\t\n", + "\t\tw_cluster = solve_intracluster_miqp(mu_cluster, Sigma_cluster, delta_k, alpha_k_val, bar_w, gamma)\n", + " \n", + "\t\tfor j, asset_idx in enumerate(cluster):\n", + "\t\t\tw[asset_idx] = w_cluster[j]\n", + "\treturn w, alpha_k\n", + "" + ], + "execution_count": 4, + "outputs": [], + "id": "cd2a5ffd" + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "language": "python", + "name": "python3" + }, + "language_info": { + "codemirror_mode": { + "name": "ipython", + "version": 3 + }, + "file_extension": ".py", + "mimetype": "text/x-python", + "name": "python", + "nbconvert_exporter": "python", + "pygments_lexer": "ipython3", + "version": "3.10.10" + } }, - { - "cell_type": "markdown", - "id": "71b7d252", - "metadata": {}, - "source": [ - "### Portfolio Optimization Problem Formulation\n", - "\n", - "The classical **mean-variance portfolio optimization** problem, introduced by Markowitz (1952), aims to minimize portfolio risk for a given level of expected return (or equivalently, balance both using a risk aversion parameter). It is formulated as:\n", - "\n", - "$$\n", - "\\min_{w \\in \\mathbb{R}^N} \\quad w^\\top \\Sigma w - \\gamma \\mu^\\top w \\quad \\text{subject to} \\quad \\mathbf{1}^\\top w = 1,\\quad 0 \\leq w_i \\leq 1 \\quad \\forall i \\in \\{1, \\ldots, N\\}\n", - "$$\n", - "\n", - "where:\n", - "\n", - "- $ w = (w_1, \\ldots, w_N)^\\top $ is the vector of portfolio weights \n", - "- $ \\Sigma \\in \\mathbb{R}^{N \\times N} $ is the **covariance matrix** of asset returns \n", - "- $ \\mu \\in \\mathbb{R}^N $ is the **expected return vector** \n", - "- $ \\gamma > 0 $ is the **risk aversion coefficient** \n", - "- $ \\mathbf{1} \\in \\mathbb{R}^N $ is a vector of ones \n", - "- The constraint $ \\mathbf{1}^\\top w = 1 $ ensures the portfolio is fully invested \n", - "- The box constraint $ 0 \\leq w_i \\leq 1 $ limits individual asset weights\n", - "\n", - "This is a **convex quadratic program** and can be efficiently solved using standard optimization methods.\n", - "\n" - ] - }, - { - "cell_type": "markdown", - "id": "a7602ae4", - "metadata": {}, - "source": [ - "This classical model is easy to solve due to its convexity, but it does not capture important practical constraints that investors often face. In real-world applications, it is common to:\n", - "\n", - "- Reduce transaction costs \n", - "- Control portfolio complexity \n", - "- Avoid excessive concentration in a few assets\n", - "\n", - "To reflect these needs, we introduce two key extensions:\n", - "\n", - "#### 1. Cardinality Constraint\n", - "\n", - "Investors may want to limit the total number of assets in the portfolio. Let $\\Delta \\in \\mathbb{N}_+$ denote the maximum allowed number of assets. This leads to the **cardinality constraint**:\n", - "\n", - "$$\n", - "\\|w\\|_0 \\leq \\Delta\n", - "$$\n", - "\n", - "Here, $\\|w\\|_0$ is the **zero-norm**, which counts the number of nonzero entries in $w$. This constraint is **non-convex** and introduces **combinatorial complexity**, requiring the use of binary decision variables.\n", - "\n", - "#### 2. Box Constraint on Asset Weights\n", - "\n", - "To prevent overexposure to individual assets, we constrain each weight to lie within a maximum threshold. Let $\\bar{w} \\in (0, 1]$ be the maximum allowable weight per asset. The constraint is written as:\n", - "\n", - "$$\n", - "w \\in [0, \\bar{w}]^N\n", - "$$\n", - "\n", - "This **box constraint** ensures no single asset can dominate the portfolio, and also helps with regulatory compliance or risk diversification.\n", - "\n", - "Together, these constraints transform the problem into a **Mixed-Integer Quadratic Program (MIQP)**, which is computationally much harder to solve than the unconstrained case.\n", - "\n", - "We can then formulate the **Cardinality-Constrained Mean-Variance (CCMV)** optimization problem as follows:\n", - "\n", - "$$\n", - "w^*(\\Delta, \\bar{w}) = \\arg\\min_{w \\in \\mathbb{R}^N} \\quad w^\\top \\Sigma w - \\gamma \\mu^\\top w \\quad \\text{subject to} \\quad \\mathbf{1}^\\top w = 1, \\quad 0 \\leq w_i \\leq \\bar{w} \\quad \\forall i, \\quad \\|w\\|_0 \\leq \\Delta\n", - "$$\n", - "\n", - "Furthermore, it can be shown that this problem is NP-hard due to the combinatorial nature of the cardinality constraint. The introduction of binary variables to enforce this constraint complicates the optimization significantly.\n" - ] - }, - { - "cell_type": "markdown", - "id": "de9a9d14", - "metadata": {}, - "source": [ - "### Solving the CCMV Problem\n", - "\n", - "The **Cardinality-Constrained Mean-Variance (CCMV)** problem is computationally challenging due to its non-convex nature and the presence of binary decisions. To address this, we adopt the **clustering-based approach** proposed in the paper.\n", - "\n", - "Specifically, we apply **hierarchical clustering** to group assets based on the similarity of their return profiles. Each cluster represents a group of highly correlated assets, and we solve the CCMV problem using the **cluster centroids** as representative assets.\n", - "\n", - "This technique offers two key advantages:\n", - "- It **reduces the dimensionality** of the problem, making it more tractable.\n", - "- It **preserves key structure** in the data by exploiting the correlation between assets within each cluster.\n", - "\n", - "By solving a reduced CCMV problem on the cluster level, we approximate the original high-dimensional problem while significantly improving computational efficiency.\n" - ] - }, - { - "cell_type": "markdown", - "id": "cbf7e6c1", - "metadata": {}, - "source": [ - "### Implementation of Hierarchical Clustering \n", - "\n", - "We now implement the hierarchical clustering approach, following the pseudo-code provided in the paper. This algorithm clusters assets based on the similarity of their return profiles by analyzing the correlation structure of the return matrix.\n", - "\n", - "The purpose of this step is to reduce the dimensionality of the portfolio optimization problem by grouping similar assets together. This allows us to simplify the later stages of solving the CCMV problem, while still preserving the core structure of the asset universe.\n", - "\n", - "The clustering is performed using a bottom-up (agglomerative) approach, with distances computed between cluster centroids in the distance space. The final partition is determined by cutting the dendrogram at the optimal level, based on the largest gap in merge distances." - ] - }, - { - "cell_type": "code", - "execution_count": 2, - "id": "95332e05", - "metadata": {}, - "outputs": [], - "source": [ - "import numpy as np\n", - "import cvxpy as cp\n", - "import matplotlib.pyplot as plt" - ] - }, - { - "cell_type": "code", - "execution_count": 3, - "id": "88e623bf", - "metadata": {}, - "outputs": [], - "source": [ - "def Hierarchical_Clustering(R: np.ndarray) -> np.ndarray:\n", - "\t\"\"\"\n", - "\tPerform hierarchical clustering on the given return matrix R.\n", - "\t\n", - "\tParameters:\n", - "\tR (np.ndarray): A return matrix.\n", - "\t\n", - "\tReturns:\n", - "\tnp.ndarray: A 1D array of clusters.\n", - "\t\"\"\"\n", - "\n", - "\t# 1. Compute the correlation amtrix tau\n", - "\ttau = np.corrcoef(R, rowvar=False)\n", - "\n", - "\t# 2. Construct the distance matrix D\n", - "\tD = 1-tau\n", - "\n", - "\t# 3. Initialize clusters C\n", - "\tC = [[i] for i in np.arange(R.shape[1])]\n", - "\n", - "\tdef d(Ci, Cj):\n", - "\t\t\"\"\"\n", - "\t\tCompute the distance between two clusters Ci and Cj.\n", - "\t\t\n", - "\t\tParameters:\n", - "\t\tCi (list): First cluster.\n", - "\t\tCj (list): Second cluster.\n", - "\t\t\n", - "\t\tReturns:\n", - "\t\tfloat: The distance between the two clusters.\n", - "\t\t\"\"\"\n", - "\t\tcent_i = 1/ len(Ci) * np.sum(D[Ci, :], axis=0)\n", - "\t\tcent_j = 1/ len(Cj) * np.sum(D[Cj, :], axis=0)\n", - "\t\tdist = (len(Ci) * len(Cj)) / (len(Ci) + len(Cj)) * np.linalg.norm(cent_i - cent_j) ** 2\n", - "\t\treturn dist\n", - "\t# 4. Repeat until only one cluster remains\n", - "\tmerge_distances = []\n", - "\tmerge_log = []\n", - "\twhile len(C) > 1:\n", - "\t\tmin_distance = float('inf')\n", - "\t\tbest_pair = (0, 0)\n", - "\t\tfor i in range(len(C)):\n", - "\t\t\tfor j in range(i+1, len(C)):\n", - "\t\t\t\tdistance = d(C[i], C[j])\n", - "\t\t\t\tif distance < min_distance:\n", - "\t\t\t\t\tmin_distance = distance\n", - "\t\t\t\t\tbest_pair = (i, j)\n", - "\t\t\n", - "\t\t# Merge the best pair of clusters\n", - "\t\ti, j = best_pair\n", - "\t\tnew_cluster = C[i] + C[j]\n", - "\t\tmerge_log.append((C[i], C[j], new_cluster, min_distance))\n", - "\t\tdel C[max(i, j)] # Remove the larger index first to avoid index issues\n", - "\t\tdel C[min(i, j)]\n", - "\t\tC.append(new_cluster)\n", - "\t\tmerge_distances.append(min_distance)\n", - "\n", - "\t# 5. Determine optimal number of clusters\n", - "\tgap = []\n", - "\tfor k in range(1, len(merge_distances)):\n", - "\t\tgap.append(merge_distances[k] - merge_distances[k-1])\n", - "\tlargest_gap = np.argmax(gap) + 1\n", - "\tthreshold = merge_distances[largest_gap]\n", - "\t# 6. Assign clusters based on the threshold\n", - "\tclusters = [{i} for i in range(R.shape[1])]\n", - "\n", - "\tfor log in merge_log:\n", - "\t\tif log[3] <= threshold:\n", - "\t\t\tCi, Cj = log[0], log[1]\n", - "\t\t\tto_merge = []\n", - "\t\t\tfor cluster in clusters:\n", - "\t\t\t\tif any(i in cluster for i in Ci) or any(i in cluster for i in Cj):\n", - "\t\t\t\t\tto_merge.append(cluster)\n", - "\n", - "\t\t\tif len(to_merge) == 2:\n", - "\t\t\t\tclusters.remove(to_merge[0])\n", - "\t\t\t\tclusters.remove(to_merge[1])\n", - "\t\t\t\tclusters.append(set(to_merge[0]).union(set(to_merge[1])))\n", - "\t\n", - "\tlabels = np.empty(R.shape[1], dtype=int)\n", - "\tfor cluster_id, asset_set in enumerate(clusters):\n", - "\t\tfor asset in asset_set:\n", - "\t\t\tlabels[asset] = cluster_id\n", - "\treturn labels" - ] - }, - { - "cell_type": "markdown", - "id": "40b512f4", - "metadata": {}, - "source": [ - "### Modeling the Clustered CCMV Problem\n", - "\n", - "Let $C = \\{C_1, \\dots, C_K\\}$ be the optimal clusters obtained from hierarchical clustering. We now introduce three key parameters:\n", - "\n", - "- Maximum allowable weight for any asset $\\bar{w}$\n", - "- Maximum number of assets per cluster $\\Delta = (\\Delta_1, \\dots \\Delta_k)$\n", - "- $\\alpha = (\\alpha_1, \\dots, \\alpha_k)$, where $\\alpha_i$ is the fraction of the portfolio allocated to cluster $C_i$\n", - "\n", - "We can then define the feasible set of portfolio weights accounting for both clustering and cardinality constraints:\n", - "$$\n", - "\t\\mathcal{A}_C(\\Delta, \\alpha, \\bar{w}) = \\left\\{ w \\in [0, \\bar{w}]^N : \\|w_{C_i}\\|_0 \\leq \\Delta_i, \\quad \\mathbf{1}^\\top w_{C_i} = \\alpha_i \\quad \\forall i = 1, \\ldots, K \\right\\}\n", - "\n", - "$$\n", - "\n", - "We can then express the clustered CCMV problem as:\n", - "\n", - "$$\n", - "\tw^*(C, \\Delta, \\alpha, \\bar{w}) = \\arg\\min_{w \\in \\mathcal{A}_C(\\Delta, \\alpha, \\bar{w})} \\quad w^\\top \\Sigma w - \\gamma \\mu^\\top w\n", - "$$\n", - "\n", - "Our goal is to determine the optimal cluster-level allocation $\\alpha_k$, the maximum number of assets $\\Delta_k$ per cluster $C_k$, and the optimal asset weights $w$ to balance allocation and enchance portfolio diversification.\n", - "\n", - "We now introduce an intra-cluster optimization problem to determine the optimal asset weights within each cluster. This is formulated as:\n", - "$$\n", - "\tw^*_{C_k} = \\arg\\min_{w \\in \\mathcal{A}_{C_k}(\\Delta_k, \\alpha_k, \\bar{w})} \\quad w^\\top \\Sigma_{C_k} w - \\gamma \\mu_{C_k}^\\top w\n", - "$$\n", - "\n", - "Where $\\Sigma_{C_k}$ and $\\mu_{C_k}$ are the covariance matrix and expected return vector for the assets in cluster $C_k$. We can also define the feasible set used above as:\n", - "$$\n", - "\\mathcal{A}_{C_k}(\\Delta_k, \\alpha_k, \\bar{w}) = \\left\\{ w \\in [0, \\bar{w}]^{|C_k|} : \\|w\\|_0 \\leq \\Delta_k, \\quad 1^Tw = \\alpha_k \\right\\}.\n", - "$$\n", - "\n", - "This means that we can write the original feasible set as:\n", - "$$\n", - "\\mathcal{A}_C(\\Delta, \\alpha, \\bar{w}) = \\bigcap_{k=1}^K \\left\\{ w : w_{C_k} \\in \\mathcal{A}_{C_k}(\\Delta_k, \\alpha_k, \\bar{w}) \\right\\}.\n", - "$$\n", - "\n", - "We now look at the objective function. To simplify things, we neglect correlations between different clusters, which allows us to write the objective function as a sum of intra-cluster objectives:\n", - "$$\n", - "f(w) \\approx \\sum_{k=1}^K f_{C_k}(w) = \\sum_{k=1}^K \\left( \\sum_{i,j \\in C_k} w_iw_j\\sigma_{ij} -\\gamma \\sum_{i \\in C_k} w_i \\mu_i \\right).\n", - "$$\n", - "\n", - "We can then approximate the optimal objective value as \n", - "$$\n", - "f^*_C(\\Delta, \\alpha) \\approx \\sum_{k=1}^K f^*_{C_k}(\\Delta_k, \\alpha_k)\n", - "$$\n", - "\n", - "where $f^*_C(\\Delta, \\alpha) = \\min_{w \\in \\mathcal{A}_C(\\Delta, \\alpha, \\bar{w})} f(w)$ and $f^*_{C_k}(\\Delta_k, \\alpha_k) = \\min_{w_{C_k} \\in \\mathcal{A}_{C_k}(\\Delta_k, \\alpha_k, \\bar{w})} f_{C_k}(w)$.\n", - "\n", - "We can then finally define\n", - "\n", - "$$\n", - "\tf^*_C(\\Delta) := \\min_{\\substack{\\Delta, \\alpha, w_i; \\\\ \\sum_{k=1}^K \\Delta_k = \\Delta; \\\\ \\sum_{k=1}^K \\alpha_k = 1}} f^*_C(\\Delta, \\alpha) = \\min_{\\substack{\\Delta, \\alpha, w_i; \\\\ \\sum_{k=1}^K \\Delta_k = \\Delta; \\\\ \\sum_{k=1}^K \\alpha_k = 1}} \\min_{w \\in \\mathcal{A}_C(\\Delta, \\alpha, \\bar{w})} f(w)\n", - "$$\n", - "\n", - "This final formulation represents a two-level optimization: the outer minimization selects how to allocate weight ($\\alpha$) and sparsity ($\\Delta$) across clusters, while the inner minimization solves for the optimal asset weights within each cluster under those constraints. Together, this balances diversification across clusters with precise selection within them.\n" - ] - }, - { - "cell_type": "markdown", - "id": "dab701ff", - "metadata": {}, - "source": [ - "### Implementation of $\\Delta$-CCMV and $\\alpha$-CCMV\n", - "\n", - "We now implement the $\\Delta$-CCMV and $\\alpha$-CCMV problems as described by the pseuducode in the paper. The first method $\\Delta$-CCMV pre-assigns the number of selected assets per cluster, while the second method $\\alpha$-CCMV pre-assigns the fraction of the portfolio allocated to each cluster." - ] - }, - { - "cell_type": "markdown", - "id": "d751606d", - "metadata": {}, - "source": [ - "#### $\\Delta$-CCMV\n", - "\n", - "We introduce a score based method to determine the number of selected assets per cluster. For each cluster $C_k$, we compute a score based on the expected return and risk of the assets within that cluster.\n", - "We assign $\\Delta_k = \\lfloor{\\Delta \\times s_k}\\rfloor$ to each cluster $C_k$ where $s_k$ is the relative importance of the cluster. In the algorithm, we compute the score as follows:\n", - "$$\n", - "s_k = \\left(\\frac{\\gamma \\mu_{C_k} - \\lambda}{2 \\sigma^2_{C_k}}\\right)_+\n", - "$$\n", - "\n", - "where $\\mu_{C_k}$ and $\\sigma^2_{C_k}$ are the expected return and variance of the assets in cluster $C_k$. $\\lambda$ is a normalization constant that ensures the scores sum to 1.\n", - "\n", - "By using scores, we can ensure that the clusters with stronger risk-adjusted preformance are more represented in the final portfolio, while still respecting the cardinality constraint." - ] - }, - { - "cell_type": "code", - "execution_count": 4, - "id": "cd2a5ffd", - "metadata": {}, - "outputs": [], - "source": [ - "from scipy.optimize import brentq\n", - "from math import floor\n", - "import cvxpy as cp\n", - "import numpy as np\n", - "\n", - "def allocate_cardinality(clusters, mu, Sigma, Delta, gamma):\n", - " \"\"\"\n", - " Allocate the number of assets (Delta_k) to each cluster based on the cluster-level \n", - " expected return and variance, as defined in Step 2 of the Δ-CCMV algorithm.\n", - "\n", - " Parameters:\n", - " clusters (list of lists): Each sublist contains indices of assets in a cluster.\n", - " mu (np.ndarray): Vector of expected returns for all assets (length N).\n", - " Sigma (np.ndarray): Covariance matrix of asset returns (shape N x N).\n", - " Delta (int): Total number of assets to include in the portfolio (cardinality constraint).\n", - " gamma (float): Risk aversion parameter.\n", - "\n", - " Returns:\n", - " list of int: A list Delta_k of length K (number of clusters), where Delta_k[i] is the \n", - " number of assets allocated to cluster i. The sum of all Delta_k equals Delta.\n", - " \"\"\"\n", - "\n", - " # Step 1: Compute cluster-level expected returns and variances\n", - " cluster_returns = []\n", - " cluster_variances = []\n", - " for cluster in clusters:\n", - " cluster_assets = np.array(cluster)\n", - " cluster_mu = mu[cluster_assets]\n", - " cluster_Sigma = Sigma[np.ix_(cluster_assets, cluster_assets)]\n", - " \n", - " expected_return = np.mean(cluster_mu)\n", - " variance = np.mean(np.diag(cluster_Sigma))\n", - " \n", - " cluster_returns.append(expected_return)\n", - " cluster_variances.append(variance)\n", - " \n", - " # Step 2: Define the score function\n", - " def score_function(lambda_, mu_k, sigma2_k, gamma):\n", - " scores = []\n", - " for mu_i, sigma2_i in zip(mu_k, sigma2_k):\n", - " raw = (gamma * mu_i - lambda_) / (2 * sigma2_i)\n", - " scores.append(max(0, raw))\n", - " return scores\n", - "\n", - " def root_objective(lambda_, mu_k, sigma2_k, gamma):\n", - " scores = score_function(lambda_, mu_k, sigma2_k, gamma)\n", - " return sum(scores) - 1\n", - "\n", - " # Step 3: Solve for lambda*\n", - " lambda_star = brentq(root_objective, -1000, 1000, args=(cluster_returns, cluster_variances, gamma))\n", - "\n", - " # Step 4: Final scores and allocation\n", - " scores = score_function(lambda_star, cluster_returns, cluster_variances, gamma)\n", - " Delta_k = [floor(Delta * s) for s in scores]\n", - "\n", - " return Delta_k\n", - "\n", - "def solve_intracluster_miqp(mu_cluster, Sigma_cluster, delta_k, alpha_k, bar_w, gamma):\n", - " \"\"\"\n", - " Solve the MIQP problem for a single cluster C_k.\n", - "\n", - " Parameters:\n", - " mu_cluster (np.ndarray): Expected returns for assets in the cluster (shape: (n,))\n", - " Sigma_cluster (np.ndarray): Covariance matrix for the cluster (shape: (n, n))\n", - " delta_k (int): Number of assets to select in this cluster.\n", - " alpha_k (float): Fraction of the total portfolio weight allocated to this cluster.\n", - " bar_w (float): Maximum allowed weight per asset.\n", - " gamma (float): Risk aversion parameter.\n", - "\n", - " Returns:\n", - " np.ndarray: Optimal weights for the assets in the cluster (length n).\n", - " \"\"\"\n", - " n = len(mu_cluster) # number of assets in the cluster\n", - "\n", - " # Step 1: Define variables\n", - " w = cp.Variable(n)\n", - " z = cp.Variable(n, boolean=True)\n", - "\n", - " # Step 2: Define the objective function\n", - " objective = cp.Minimize(cp.quad_form(w, Sigma_cluster) - gamma * mu_cluster @ w)\n", - "\n", - " # Step 3: Define the constraints\n", - " constraints = [\n", - " cp.sum(w) == alpha_k,\n", - " w >= 0,\n", - " w <= bar_w * z,\n", - " cp.sum(z) <= delta_k\n", - " ]\n", - "\n", - " # Step 4: Solve the problem\n", - " problem = cp.Problem(objective, constraints)\n", - " problem.solve(solver=cp.GUROBI if cp.GUROBI in cp.installed_solvers() else cp.ECOS_BB)\n", - "\n", - " # Step 5: Return the solution\n", - " return w.value\n", - "\t\n", - "\n", - "def delta_ccmv(R, mu, Sigma, Delta, bar_w, gamma):\n", - "\t\"\"\"\n", - "\tSolve the CCMV problem by pre-assigning the number of assets in each cluster.\n", - "\t\n", - "\tParameters:\n", - "\tR (np.ndarray): A return matrix.\n", - "\tmu (np.ndarray): Expected returns.\n", - "\tSigma (np.ndarray): Covariance matrix of returns.\n", - "\tDelta (int): Total number of assets to include in the portfolio (cardinality constraint)\n", - "\tbar_w (float): Maximum allowed weight per asset\n", - "\tgamma (float): Risk aversion parameter.\n", - "\t\n", - "\tReturns:\n", - "\tw: Optimal portfolio weights\n", - "\talpha: cluster allocation proportions\n", - "\t\"\"\"\n", - "\n", - "\t# Step 1: Collect the clusters from R\n", - "\tclusters = Hierarchical_Clustering(R)\n", - "\n", - "\t# Step 2: Delta pre-assignment\n", - "\t\n", - "\tDelta_k = allocate_cardinality(clusters, Delta)\n", - " \n", - "\t# Step 3: Optimization\n", - " \n", - "\talpha_k = [1 / len(clusters) for _ in clusters]\n", - "\tw = np.zeros(R.shape[1])\n", - "\tfor i, cluster in enumerate(clusters):\n", - "\t\tmu_cluster = mu[cluster]\n", - "\t\tSigma_cluster = Sigma[np.ix_(cluster, cluster)]\n", - "\t\tdelta_k = Delta_k[i]\n", - "\t\talpha_k_val = alpha_k[i]\n", - "\t\t\n", - "\t\tw_cluster = solve_intracluster_miqp(mu_cluster, Sigma_cluster, delta_k, alpha_k_val, bar_w, gamma)\n", - " \n", - "\t\tfor j, asset_idx in enumerate(cluster):\n", - "\t\t\tw[asset_idx] = w_cluster[j]\n", - "\treturn w, alpha_k\n" - ] - } - ], - "metadata": { - "kernelspec": { - "display_name": "Python 3", - "language": "python", - "name": "python3" - }, - "language_info": { - "codemirror_mode": { - "name": "ipython", - "version": 3 - }, - "file_extension": ".py", - "mimetype": "text/x-python", - "name": "python", - "nbconvert_exporter": "python", - "pygments_lexer": "ipython3", - "version": "3.10.10" - } - }, - "nbformat": 4, - "nbformat_minor": 5 -} + "nbformat": 4, + "nbformat_minor": 5 +} \ No newline at end of file diff --git a/README.md b/README.md index b0b0a3f..6d6149a 100644 --- a/README.md +++ b/README.md @@ -30,7 +30,7 @@ $$ - ✅ Solves the **cardinality-constrained minimum-variance** problem - ✅ Compares with the unconstrained Markowitz solution -- ✅ Uses `mip` for mixed-integer quadratic programming +- ✅ Uses SciPy to solve small cardinality-constrained examples by enumerating selected asset subsets - ✅ Visualizes how limiting asset count affects portfolio composition and performance - ✅ Prepares groundwork for future QUBO or quantum annealing implementations @@ -48,10 +48,11 @@ $$ - Python 3.x - NumPy +- SciPy - Matplotlib -- `mip` (Mixed Integer Programming) Install dependencies: ```bash -pip install numpy matplotlib mip +pip install -r requirements.txt +``` diff --git a/requirements.txt b/requirements.txt new file mode 100644 index 0000000..f23634b --- /dev/null +++ b/requirements.txt @@ -0,0 +1,3 @@ +numpy +scipy +matplotlib