Analytical Development of the Efficient Frontier of Portfolio and Electricity Generation in Mexico

One of the programs for electricity generation forecast by the Centro Nacional de Energía (CENACE-National Center for Energy Control) based on cost structure by technology type and existing capacity is constituted by the optimization model based on the branch and bound algorithm (Branch and Bound). The PLEXOS option of this algorithm is solved by linear programming (LP), that is to say, without the existence of integrity conditions over the variables [1].

Taking into account installed capacity, in 2017 structure according to electricity generation technology behaved as follows: for combined-cycle (CC) 50%; for thermoelectric power (TCC) 19.4%; for carboelectric power (CAR) 9.3%; of hydroelectric power (HIDRO) 9.7%; of nuclear generation (NUC) 3.3%; of wind power (EOLO) 3.2% and, of geothermal and biofuel technology 4.9%. For 2032, CENACE's energy mix is estimated as follows: of combined cycle 41.9%; of conventional thermoelectric power 10.1%; of carboelectric power 3.2%; of hydroelectricity 11.4%; geothermal energy 1.5%; of nuclear power 4.4%; of wind power 14.6%, and the remaining being a combination of geothermal, solar photovoltaic, biomass and efficient cogeneration.

Purpose of this work is to contrast cost structure per forecasted technology for the year 2032 using the least-cost method (currently used by CENACE) with option of mean-variance approach for Mexico. To this effect, the efficient frontier of asset portfolio will be generated according to the return-risk pair. This description is due to the approach goal which is based on obtaining the efficient technology portfolio implying to expose society to the risk minimum level necessary related to electricity generation.

In finance, investors use portfolio theory to minimize the risk and maximize the return of their portfolio by diversification. Some theories can be used in electricity generation planning to evaluate the addition of each new alternative on the basis of its contribution to the overall risk and cost of the generation mix and not on its stand – alone cost [5].

There are electricity generation technologies showing high correlations, according to their fuel and/or operation and maintenance costs among them (more substitutes) such as CC gas plants as well as those of fossil fuel and coal; likewise, other technologies with almost nil correlations (more complementary) can be identified, such as those using fossil fuels with those operating based on nuclear energy or with solar and wind technologies. Besides, there is also a possibility of technologies with a negative correlation, such as it is the case between Biomass and fuel derivatives [12].

These considerations are the basis for implementing optimization processes through theories involving mean-variance use. Having this information, one of the most interesting questions in the electricity sector area is the following: what changes may be feasible in the energy mix of a country or region allowing to win in one of both dimensions (risk-return) of optimal portfolio without harming the other? Of course, among all efficient portfolios, the optimal choice depends on preferences of the energy policy consultant [17].

3.1 Set of efficient portfolios when all values are at risk

On the Markowitz's portfolios optimization theory formalism, the main elements are a set of m assets. From each one of them arises a series in time of T data (measured in investment units); in the case of investment portfolios, they correspond to assets' closing prices. Thus, from which we obtain returns T-1, following Eq. (1):

$r_t=\left(P_t-P_{t-1}\right) / P_{t-1}($ dimensionless $)$               (1)

These returns T-1 form a series in time having those same elements number.

3.2 Portfolios optimization problem

Markowitz's attainment was to relate standard deviation with the risk, giving to the latter the task of assigning a quantitative measure. The problem of reducing the risk of a portfolio, given an expected return value, can be mathematically written as an optimization problem with m-dimensional constraints that is given by:

$\min _{\underline{x}_p} F\left(\underline{x}_p\right)$               (2)

s.t. $E\left(\underline{x}_p\right)=E_p$               (3)

$\sum_{i=1}^m x_i=1$               (4)

Solution to this problem is $\underline{x}_p^*$, in such a way that it complies with being the minimum of all values of objective function F and additionally with the constraints that $E\left(\underline{x}_p^*\right)=E_p$ and $\sum_{i=1}^m x_i^*=1$.

We can change the optimization problem with m - dimensional constraints to a problem with no constraints m + 2 - dimensional, using the Lagrange multipliers:

$\min _{\left(\underline{x}_p, \lambda_1, \lambda_2\right)} L:=F\left(\underline{x}_p\right)-\gamma_1 g_1\left(\underline{x}_p\right)-\gamma_2 g_2\left(\underline{x}_p\right)$               (5)

where, g₁ and g₂ are defined using problem constrains:

$g_1\left(\underline{x}_p\right):=E\left(\underline{x}_p\right)-E_p$               (6)

$g_2\left(\underline{x}_p\right):=\sum_{i=1}^m\left(\underline{x}_p\right)_i-1$               (7)

3.3 Solving the problem with no constrains

Conditions necessary to find a stationary point for the problem with no constrains are given by:

$\nabla_{\left(x, \lambda_1, \lambda_2\right)} L=0$               (8)

The $\nabla_{\left(\underline{x}, \lambda_1, \lambda_2\right)}$ operator means that we will perform m+2 partial derivatives of L, we will accommodate them as components of a column vector and we equalize them to a column vector with zeros as components:

$\frac{\partial}{\partial x_1} L=0 \cdots \frac{\partial}{\partial x_m} L=0 \cdot \frac{\partial}{\partial \lambda_1} L=0, \frac{\partial}{\partial \lambda_2} L=0$               (9)

To find the stationary point $\left(\underline{x}_p^*, \lambda_1^*, \lambda_2^*\right)$ is necessary to solve the algebraic equations system resulting from the solution of previous partial derivatives, which can be explicitly written as follows:

$0=\sum_{i=1}^m x_j \sigma_{i j}-\gamma_1 E_i-\gamma_2 \cdot 1, i=1, \ldots m$               (10)

$0=E-\sum_{i=1}^m x_i E_i$               (11)

$0=1-\sum_{i=1}^m x_i$               (12)

Please note that we have omitted bars and asterisks from our previous notation so that $\left(\underline{x}_p^*\right)_j=x_j, \gamma_1^*=\gamma_1 y \gamma_2^*=\gamma_2$.

Here we can see that there are m+2 equations with m+2 unknowns instead of a system of partial differential equations. For convenience we will go to the matrix form of the problem:

$\left(\sigma \underline{x}^T\right)_i-\left(\gamma_1 \underline{E}\right)_i-\left(\gamma_2 \underline{1}\right)_i=0, i=1, \ldots, m$               (13)

$0=E-\underline{x}^T \underline{E}$               (14)

$0=1-\underline{1}^T \underline{x}$               (15)

The weights vector (column) $\underline{x}=\left(x_1, x_2, \ldots, x_m\right)^T$ is the vector that contains the percentages (which will be a main part of the optimization problem together with the Lagrange multipliers $\lambda_1$ and $\lambda_2$. We also specify that $\underline{x}^T=\left(x_1, x_2, \ldots, x_m\right)$ is a row or line vector of m inputs with the same components as $\underline{x}$. Here $\underline{E}=\left(\mu_1, \ldots, \mu_m\right)^T$ is the column vector containing the averages or means of historical data. On the other hand, $\underline{1}=(1, \ldots, 1)^T$ is the column vector of m components that only has entries with value 1.

Making use of matrix notation to find the solution to the equations system:

$\sigma \underline{x}-\gamma_1 E-\gamma_2 \underline{1}=\underline{0}$               (16)

where, $\underline{0}^T=(0, \ldots, 0)^T$ is a column vector with zeros in all its inputs (components). We define the inverse matrix of $\sigma$ as the one having the v_ij components:

$\left[\sigma^{-1}\right]_{i j}=v_{i j}$               (17)

Analytically speaking (meaning, beyond computers field) for the inverse matrix of σ, to exist, sufficient condition that must be met is that determinant must be different to zero.

Therefore, by multiplying Eq. (16) by the inverse matrix we have:

$\sigma^{-1} \sigma \underline{x}-\gamma_1 \sigma^{-1} \underline{E}-\gamma_2 \sigma^{-1} \underline{1}=\sigma^{-1} \underline{0}$               (18)

note that $\sigma^{-1} \sigma=1_{m \times m}$, where $1_{m \times m}$ is the identity matrix that has ones on the diagonal and zeros in all the remaining entries. By multiplying inverse matrix by a vector, we obtain another same dimension vector as the initial one. And finally, we observe the simple multiplication $\sigma^{-1} \underline{0}=\underline{0}$ on the equation's right side.

Finally, the result of the above equation is:

$\underline{x}-\gamma_1 \sigma^{-1} \underline{E}-\gamma_2 \sigma^{-1} \underline{1}=\underline{0}$               (19)

note that if $\gamma_1$ and $\gamma_2$ were known, then we would immediately have the solution values since vector $\underline{E}$ and matrix $\sigma^{-1}$ are known or calculable from the problem data. Therefore, our task from now on is to find $\gamma_1$ and $\gamma_2$.

We can write for the k-th component of Eq. (19) as:

$x_k=\gamma_1 \sum_{j=1}^m v_{k j} E_j+\gamma_2 \sum_{j=1}^m v_{k j}, k=1, \ldots m$               (20)

Multiplying the previous equation by E_k and adding on k we have:

$\sum_1^m x_k E_k=\gamma_1 \sum_1^m \sum_1^m v_{k j} E_j E_k+\gamma_2 \sum_1^m \sum_1^m v_{k j} E_k$               (21)

which in matrix notation is:

$\underline{x}^T \underline{E}=\gamma_1 \underline{E}^T \sigma^{-1} \underline{E}+\gamma_2 \underline{1}^T \sigma^{-1} \underline{E}$               (22)

on this equation's right side we will use the definition of the return expected value of a portfolio in its matrix form, $\underline{x}^T \underline{E}=E_p$, where E_p is a problem's data, thus:

$E_p=\gamma_1 A+\gamma_2 B$               (23)

where two of the three Merton scalars have been defined:

$A=\underline{E}^T \sigma^{-1} \underline{E}, B=\underline{1}^T \sigma^{-1} \underline{E}$               (24)

which can be calculated from the variance-covariance inverse matrix and from the means vector; therefore, they are values that can be considered problem's data. On the other hand, if we carry out the contraction operation, that is the sum of all the vector components x_k, we obtain:

$\sum_1^m x_k=\gamma_1 \sum_1^m \sum_1^m v_{k j} E_j+\gamma_2 \sum_1^m \sum_1^m v_{k j} E_k$               (25)

This also has its matrix form, represented by:

$\underline{1}^T \underline{x}=\gamma_1 B+\gamma_2 C$               (26)

where we can see that the third scalar, obtained from Merton is:

$C \equiv \sum_1^m \sum_1^m v_{k j}$               (27)

Left side of this Eq. (26) equals 1 as it is the budget constraint to which the problem is subject, and therefore:

$1=\gamma_1 B+\gamma_2 C$               (28)

Thus, a system of two equations with two unknowns is formed $\left(\gamma_1, \gamma_2\right)$ whose solution is:

$\gamma_1=\frac{(E C-B)}{D}, \gamma_2=\frac{(A-E B)}{D}$               (29)

This solution can be seen immediately if we use Cramer's formula and if, besides, for this solution to exist, we ask that the discriminant $B C-A^2=D \neq 0$.

Once obtained $\left(\gamma_1, \gamma_2\right)$ we substitute in 1 system to obtain the solution:

$\begin{aligned} & x_k=\frac{E \sum_1^m v_{k j}\left(C E_j-A\right)+\sum_1^m v_{k j}\left(B-A E_j\right)}{D} \\ & k=1, \ldots, m \\ & \end{aligned}$               (30)

3.4 Efficient frontier

Graphically, the efficient frontier is described by all those points on plane $\sigma-E$ where for a certain given value of E_p we find minimum $\sigma_p$ that a linear combination with weights $\underline{x}_p$ can give.

To obtain the geometric place in R² with axis of the abscissae equal to the risk $\sigma$ and the ordinates axis equal to portfolio return E_p, we construct an expression that relates these two quantities multiplying (10) by x_i and adding of $i=1, \ldots, m$:

$\sum_{i=1}^m \sum_{j=1}^m x_i x_j \sigma_{i j}=\gamma_1 \sum_{\substack{1 \\ i=1,...m}}^m x_i E_i+\gamma_2 \sum_1^m x_i$               (31)

Left side of this expression is the portfolio variance definition, that is to say, $\sigma^2\left(x_p\right)$, while on the right side we obtain $\gamma_1$ multiplying by E_p which is the expected value of the portfolio return and $\gamma_2$, the budgets’ equation which is 1. From which we obtain:

$\sigma^2\left(x_p\right)-\frac{1}{C}=\frac{C}{D}\left(E_p-\frac{A}{C}\right)^2$               (32)

If C and D are greater than zero, then the right side of Eq. (32) will always be positive since it is written in terms of $\left(E_p-\frac{A}{C}\right)^2 \geq 0$. When $E_p=\frac{A}{C}$ this term will be zero and therefore $\sigma^2\left(x_p\right)=\frac{1}{C}$.

We define $\underline{E} \equiv \frac{A}{C}$ as the minimum expected value of a portfolio that has the minimum variance $\underline{\sigma}^2 \equiv \frac{1}{C}$.

Defining $\underline{x}_k$ so that it is the proportion of the minimum variance portfolio invested in the asset $k^{t h}$, then it is necessary to redefine:

$\underline{x}_k=\frac{\sum_1^m v_{k j}}{C}, k=1, \ldots, m$               (33)

However, it is usual to present the frontier in the standard deviation plane of the mean instead of the mean-variance plane, so its geometric representation is represented will remain as:

$\frac{\sigma^2}{\frac{1}{C}}-\frac{\left(E-\frac{A}{C}\right)^2}{\frac{D}{C}}=1$               (34)

According to analytical geometry, this equation represents a hyperbola that meets and confirms conditions of an efficient portfolio with minimum variance in $E \equiv \frac{A}{C}$ and a minimum variance of $\underline{\sigma}^2 \equiv \frac{1}{C}$.

3.5 Non-negative weights

According to literature [18], both the optimal theory model and the capital assets pricing model [19] indicate that portfolios generated with positive investment percentages will have a characteristic of being efficient within the mean-variance relationship.

However, these types of portfolios (all the positive weights) are not so common to obtain, especially if data from assets' historical series are used. In the best-case scenario, their results indicate that, only it is possible to identify a single segment of the efficient frontier where all are no-negative weights.

Written parametrically, Markowitz can be transcribed as follows:

$\max \left\{t E^{\prime} x-\frac{1}{2} x^{\prime} \sigma x \mid 1^{\prime} x=1\right\}$               (35)

where, t is a scalar parameter, x is the portfolio's n-vector of weights, and $1^{\prime} x=1$ is the budget constraint. Solution to this problem in terms of parameter t for x and t is:

$x(t)=\frac{\sigma^{-1} 1}{C}+t\left[\sigma^{-1}\left(E_p-\frac{1 A}{C}\right)\right]$               (36)

$\lambda(t)=-\frac{1}{C}+\frac{t A}{C}$               (37)

where, A, B and C are the Merton scalars from the previous section. Now let us define:

$h_0=\frac{\sigma^{-1} 1}{C}$               (38)

$h_1=\sigma^{-1}\left(E_p-\frac{1 A}{C}\right)$               (39)

So, portfolio weights can be expressed as:

$x(t)=h_0+t h_1$               (40)

It is clearly seen from this equation if we observe each component separately, then:

$x_i(t)=\left(h_0\right)_i+t\left(h_1\right)_i, i=1, \ldots, m$               (41)

that if $h 1_i>0$, then $x_i(t)$ increases as t raises and it will be positive since $t \geq\left(h_0\right)_i /\left(h_1\right)_i$.

It is possible to determine a maximum dimension and a minimum dimension for t's values in such a way that we can guarantee existence of positive weights:

$t_l=\max \left\{-\left(h_0\right)_i /\left(h_1\right)_i \mid\right.$ for all $i$ with $\left.\left(h_1\right)_i>0\right\}$               (42)

$t_u=\min \left\{-\left(h_0\right)_i /\left(h_1\right)_i \mid\right.$ for all $i$ with $\left.\left(h_1\right)_i<0\right\}$               (43)

From where it is not difficult to conclude that:

$x_i(t) \geq 0$ for all $i$ with $\left(h_1\right)_i>0$ and for all $t \geq t_l$               (44)

$x_i(t) \geq 0$ for all $i$ with $\left(h_1\right)_i<0$ and for all $t<t_\mu$               (45)

note that for $x(t)$ to be positive in all its components for some t it is necessary that:

$\left(h_0\right)_i \geq 0$ for all $i$ so that $\left(h_1\right)_i=0$               (46)

This implies that $x(t)>0$ alone and if $t_l \leq t \leq t_u$ and the above constraint is complied with literature [18]. From here it is easy to obtain a strong conclusion that if $t_l>t_u$ there are no positive components of $x(t)$.

3.6 Used information

Total levelized annual cost includes three areas: a) cost per investment; b) cost per fuel and c) cost per operation and maintenance (includes cost for water consumption). From 1992 to 2011, information was provided by Taboada González et al. [20]. This information is in Mexican pesos and it was developed with information from COPAR (Spanish acronym of Costs and Reference Parameters for Investment Project Formulation in the Electrical Sector). From 2012 to 2016, the same Taboada procedure was followed except that information was stated in dollars, so the peso-dollar exchange rate was applied to convert to pesos the annual averages obtained from each generation technology. The calculation is based on the COPAR's 2012 and 2015 documents.

Going deeper into the proposed analysis, it can be observed that for the year 2017 and taking into account installed capacity, structure per electricity generation technology behaved as follows: 50.20% for combined-cycle (CC); 19.4% for thermoelectric power (TCC); 9.3% for carboelectric power (CAR); 9.7% of hydroelectric power (HIDRO); 3.3% of nuclear generation (NUC); 3.2% of wind power (EOLO) and; 4.9% of geothermal technology and biofuel. For this combination, return on portfolio was 10.8% with an associated risk of 12.56%.

In this same vein, for the year 2032, CENACE's energy mix is forecasted as follows: 41.9% of combined cycle; 10.1% of conventional thermoelectric power; 3.2% of carboelectric power; 11.4% of hydroelectricity; 1.56% geothermal; 4.4% nuclear; 14.6% of wind power and, the remainder, a combination of geothermal, solar photovoltaic, biomass and efficient cogeneration. With this combination, return on portfolio rises to 14.9% but with an associated risk of 12.59%. It is noted that for a same level of risk associated to the energy mix, profitability maximization (or cost decrease) is lower.

However, for scenario C, it is possible to identify that the proposed energy combination is feasible in relation to the future projections made by CENACE, the risk-return levels are acceptable, and it also weighs a 250% increase in clean energy generation compared to scenarios A and B. The aforementioned would help preserve Mexico's energy sovereignty while striving to reduce pollution levels in the environment.