A Long-Term Evaluation on Transmission Line Expansion Planning with Multistage Stochastic Programming

The purpose of this paper is to apply multistage stochastic programming to the transmission line expansion planning problem, especially when uncertain demand scenarios exist. Since the problem of transmission line expansion planning requires an intensive computational load, dual decomposition is used to decompose the problem into smaller problems. Following this, progressive hedging and proximal bundle methods are used to restore the decomposed solutions to the original problems. Mixed-integer linear programming is involved in the problem to decide where new transmission lines should be constructed or reinforced. However, integer variables in multistage stochastic programming (MSSP) are intractable since integer variables are not restored. Therefore, the branch-and-bound algorithm is applied to multistage stochastic programming methods to force convergence of integer variables.In addition, this paper suggests combining progressive hedging and dual decomposition in stochastic integer programming by sharing penalty parameters. The simulation results tested on the IEEE 30-bus system verify that our combined model sped up the computation and achieved higher accuracy by achieving the minimised cost.


Introduction
Development of renewable energy has been expedited by the global effort to reduce greenhouse gas emissions. These include emissions from fossil fuels used in transportation, which are declining as electric vehicle use increases [1][2][3]. The global move to electric vehicles and new facilities to support renewable resources now requires much higher electricity consumption and generation capability. If transmission line capacity fails to catch up to the growth of them, the failure will cascade and cause shortages in the grid requiring expensive repairs for the power system network [4].
Transmission line expansion planning (TLEP) has been suggested to avoid future shortages. Economic effects of probable transmission line failures are compared against the investment budget spent on new transmission lines [5]. The TLEP is usually implemented by system planners to analyse old transmission lines that would be uneconomical in a long-term evaluation of the growth of electricity demand [6]. In the evaluation process, the planner decides locations for a new construction from the present point of view to optimise the future cost. Thus, optimal decisions for the planner are conveyed by considering probable uncertain scenarios of the growth of demand and representing the most suitable candidate among all possible realisations. To find optimal expansions in the TLEP process, the planner reconciles system reliability and energy economy over a long-term time horizon. The typical optimisation formulation of TLEP determines a minimised cost for organising transmission lines and operating generators. In this regard, the investment budget for expansions can be found by optimising the use of resources. Considering more detailed scenarios of future demand usually helps planners to mitigate this uncertainty and yields lower costs than only considering the worst-case [7].
Many stochastic optimisation processes such as the robust method, the chance constraint method and multistage stochastic programming (MSSP) have been developed to account for uncertainty in the power system [8][9][10]. Typically, they predict the impact of future events by calculating the value of their possibility as costs of the objective function and compare them to the current decisions. MSSP shows an outstanding performance among the three methods but usually requires many scenarios and variables depending on time horizons [11]. Considering many scenarios and variables helps to verify results, but the size of the computational load becomes intractable. Therefore, decomposition methods are used to break the original problem into several smaller problems, called subproblems. However, the decomposition is followed by the coupling constraint added between separated subproblems. The constraints are treated as a price for solving them independently so that one subproblem can consider the others' decisions as they should pay for not having a unified value [12].
Meanwhile, there are two types of decomposition methods to consider for uncertain electricity demands in MSSP. One is the primal decomposition, represented by Bender's method and the L-shape method [13], which divides the problem into two parts: the master and the slave problem. Yet, it is hard to calculate more than three stages, since too many master and slave pairs are generated as the time span is extended, so it is not appropriate for a long-term evaluation process. The second, dual decomposition (DD), horizontally separates the scenario tree by its time horizon, so that it is inherently possible to consider a long-term horizon since subproblems are sufficiently smaller than the original. Instead, subproblems for the DD method could generate inconsistent results of subproblems since they also separate decision variables, thus the coupling constraint is added to obtain the unified solution. The coupling constraint is defined as a nonanticipativity constraint (NC) in DD methods, which explains the indistinguishable state of variables in which all separations should be converged to the one nodal solution. However, optimality of expansions is not achieved due to the independence of subproblems when integer variables are decomposed. Indeed, some decision variables used to construct a new line in the TLEP problem are comprised by integers to define a disjunctive decision; therefore, the system planner cannot decide on the construction when results diverge.
To solve the variable convergence issue, various DD methods have been developed and a few of them are adopted in this paper. Progressive hedging (PH) can transform the NC into a Lagrangian function then heuristically approach the optimal based on the proximal operator and augmented Lagrangian [14,15]. PH can quickly obtain a reliable solution, but its convergence is not guaranteed when integer variables are involved, and the integrated solution could be infeasible if the convergence rate is not sufficient. A generalised process for the dual decomposition for stochastic integer programming (DDSIP) in [16] was developed to consider integer variables, and it solves the Lagrangian dual problem of the NC [16]. Dual variables for the NC are obtained by the proximal bundle (PB) method [17], which determines a delegate among candidates from subproblems only if it is plausible. However, it is slow because of the branch-and-bound (BB) process, which bounds constraints for integer variables and usually takes a long time to calculate branched cases.
Recent studies have consistently reported improved DD methods, and numerical results have proved reliable. Studies in [12,18] mitigated the penalty of the NC by formulating it in various ways, which verified that a better matrix formulation of the NC can improve the result. Studies in [19,20] suggested methods to transform dual functions. Since the NC can be relaxed as a dual function, better objective bounds could be obtained in those studies by reversing the approach to the formulation to solve dual variables. Studies in [21][22][23] developed methods to support the heuristic approach inherent in DD methods to achieve high convergence speed and accuracy. However, those improvements were not properly implemented in practical studies of power system. The TLEP in [24], where a long-term scenario for transmission lines was decomposed and bundled according to its similarity, also in [8,25], where the DD was introduced to make generation resource schedules with many integer variables, could have been assessed better if they adjusted their DD methods according to upper and lower bounds of results and the consistency of the results was verified over a long-term time span if methods were generalised regardless of the program size. Numerical results reported in [26,27] achieved the cost-effectiveness of plans through the methodological improvement, enhancing the convergence of stochastic variables.
In this paper, a long-term evaluation problem by MSSP is proposed for the TLEP. Typical TLEP formulations are used, where decisions on new construction and reinforcement of existing transmission lines can be made by considering DC optimal power flow (OPF). However, we extend the problem in a form of MSSP to involve 30 years (with six stages) in the plan. As the time stage increases, binary variables devastate the solving process to obtain unified results and slow down the process. By adjusting DD methods in various cases studies with different sizes, we verify that stochastic variables can deliver varying convergence and expansion results with respect to algorithm. With the analysis on those spectral outputs, general forms of stochastic optimisation of DD methods are examined, and the solving process is unfolded to construct a generalised optimisation environment on varying sizes. We combine the PH and PB methods where the penalty value in PH is concatenated as a warm-start value of the PB method since they have an equal target value. Exchanging the penalty value could speed up the medium-iteration convergence and improve the quality of the solution. This contribution can make a combined method less affected by the program size. The expansion results are compared in three aspects: simulation time, nonanticipativity and objective cost. The nonanticipativity results are delivered to rectify and classify the failed results. Thus, the method termination and the feasibility are separately considered in a scope on the results. Moreover, with modified IEEE 30-bus system, the overall optimisation formulation and expansion results are given over a 30-year time span. We observed that a methodological improvement enables cost-effective evaluation, finding the best with the coherency among the feasible candidates obtained by improved boundaries of objective functions for the NC.
The organisation of this paper is summarised as follows. In Section 3, the overall process for solving MSSP with mixed-integer variables is introduced. In Section 2, TLEP problems are modelled and are assessed from the viewpoints of economy and reliability. In Section 4, the numerical results of the TLEP are presented and analysed. In Section 5, conclusions are drawn.

Transmission Line Expansion Planning Modelling
In this subsection, we introduce our multistaged transmission line expansion planning formulation.

Objective Function of TLEP
The formulation represents one subproblem, where stage is allocated to a five-year investment period. Symbols used to formulate the TLEP problem and their descriptions are presented in Table 1. The objective function of TLEP is structured as follows: In (1), we take into account the following terms in our objective functions as a view of power system planner: transmission line expansion cost (investment cost), power generation cost (operation cost) and the cost of probable load curtailments due to single line outage (load shedding cost). The investment cost is further divided into new construction cost and reinforcement cost of existing transmission lines. We only consider peak load scenarios. K c/r is a unit costs' vector of line investments, which is the hourly long-run marginal cost (LRMC) of transmission line capacity. Its dimension is $/MWh, which balances the finance scale with operation and investment cost. The investment cost is calculated by dividing annuitized LRMC by the number of hours in a year. For this parameter, discount rates can be utilised [28], but in this paper we do not consider the present value of annual cost. K g/s is the unit cost of generation and outage, respectively, whose dimension is also $/MWh.

Constraints of TLEP
For the constraints of TLEP problems, first, a balance between power supply and demand should be obtained. This constraint can be represented as follows: In (2), the sum generator outputs PG and load shedding amount PS should be equal to the sum of power flow f whose direction is out of each bus, and power load PD. Second, constraints on expansion decision u are enumerated as follows: In (3)-(5), the DC-OPF constraints are considered. Power flows f are defined by the phase angles of the buses θ and the reactance of the transmission lines X. In (4) and (5), the power flow of newly constructed lines is constrained by a disjunctive decision u with a large number Q. In this constraint, power flow of newly constructed lines and expansion decisions are linked.
In (6) and (7), power flows in all existing lines and expansion candidates abide by the maximum capacity of line loading. This capacity can be increased depending on expansion decisions.
In (8) and (9), construction decisions are expressed as binary variables and reinforcement decisions are expressed as integer variables. Finally, in (10) and (11), the sum of expansion is constrained by one, meaning transmission lines can only undergo one construction event. Constraints on the generator output PG and load shedding amount PS are enumerated as: In (12), the maximum capacity and minimum generation output is considered to decide generation output. In (13), the load shedding amount of each bus should be within the load of each bus. Finally, in (14), the generator's cost function is simplified into piece-wise linear functions for constant a and b with H pieces. Replaced transmission lines include old lines whose regular lifetimes or life expectancies are less than their initial operation periods plus the panning time horizon.

Decomposition Method Formulations
In this section, we introduce two dual decomposition (DD) methods, PH and DDSIP, and their theoretical backgrounds. PH and PB methods propose a pathway to the optimal solution, recursively updating the decisions. The BB process enforces the convergence of integer variables.

Basic Form of MSSP and Nonanticipativity Constraint
The most intuitive and intrinsic approach to represent uncertainty in possible realisations is to formulate a finite number of scenarios with relative probabilities. This is the basic form of the MSSP called the extensive-form (EF) problem and requires the scenario tree which is a set for bundles of scenarios. A scenario tree is structured by branching to possible realisations from a precedent scenario, which gives conditional probability for its succeeding scenario, as depicted in Figure 1 (left). The basic form of the MSSP problem for the scenario tree that structured over the time horizon T stages with decision variables x = [x 1 , x 2 , ..., x T ] T can be defined as follows: In (15), a simplified form of MSSP problem for TLEP is given, where f (x) is the objective function, and G is a simple expression of a set of all constraints involved in the problem for x t . The operator E[·] represents expectation of costs over related scenarios, and the scenario tree Ξ consists of all nodes of scenarios ξ t , such that ξ t ∈ Ξ. The objective function and constraints in MSSP can be decomposed into subproblems. Therefore, f (x s ) and G s represent the objective function in (1) and constraints (2)- (14), which can be extended to the multistage problem by the time span T.
In (16), the representative formulation of DD is presented, where x s is set of variables for each scenario subproblem, and ξ s is the corresponding scenario, and separated variables x s is equal to the node solutionx, which is the NC. The structure of scenario tree and subproblems are presented in Figure 1, where DD can divide the scenario tree horizontally so the number of subproblems is determined by the number of branches. In the same figure, if the NC forx b , which is the solution of node b is satisfied in the corresponding subproblems, then it is also satisfied for the second stage variable in subproblems s 1 and s 2 , such that x s 1 ,2nd = x s 2 ,2nd =x b . However, NC cannot be presupposed since x s is the optimal solution of MSSP. DD methods add a relaxed form of NC in the objective function in (16) to mitigate the constraint.
The Lagrangian function in (17) can be solved for PH and PB methods by using heuristic approaches as follows: 1.
The PH designates implementablex by the projection but derives λ S proximally. 2.
The PB solves the dual problem of (17) to get a representative for λ s and verify its degree of improvement.

Progressive Hedging Method
We provide formulations of the PH from the Lagrangian function. Independent scenario projection to the nonanticipativity solution space is assumed in the PH. The projection operator averages the subproblem solutions by considering probabilities of scenarios. The projected solution can be regarded as optimal when the NC is satisfied for all nodes in the scenario tree. An implementable solutionx can be generated by the projection.
The assumption in (18) that the optimal is the average of subproblem solutions inherently makes the NC solution space N . The solutionx ∈ N deduces the direction toward the optimal solution, and it is aligned by the probabilities ρ s . Hence, the alignment to N is possible just once so that Proj 2 = Proj and then Equation (19) is satisfied.
Decision variables in PH are updated as follows. First, in the middle of iterative process of the PH, subproblems untied from the scenario tree are solved and saved in the subproblem space S. Second, node solutions are updated by the quadratic formulation of the augmented Lagrangian. Each subproblem in (17) is transformed to be quadratic since a penalty constant ω of the NC is updated by the proximal operator for subdifferential of the dual function, ∂ λ D(λ) = (x s −x), with step size γ [15].
To update the dual variable λ s , the PH uses constant penalties for the NC in (20), which represents the distance between S and N such that ω s = x s −x. The dual variables satisfy the NC in (21) because of (19), only if the initial value for penalties is zero, and make an even boundary for the dual variables as shown in Figure 2. In short, ω represents the NC penalties according to the projected solutions, and it is orthogonal to the Proj as in (21) such that S ⊥ N .
Finally, the augmented Lagrangian is derived from iterative updates of penalties in (20) and the Lagrangian function in (17).
The first order necessity condition to find the optimal solution for (22) yields the augmented Lagrangian as follows: The overall process is summarised in Algorithm 1.

Algorithm 1 Progressive Hedging Algorithm
Convergence strategy comparison between progressive hedging (left) and proximal bundle (right).

Proximal Bundle Method
DDSIP uses PB method to solve subproblems and bounds integer variables after PB method is terminated. The PB method solves the dual problem of (17) by sharing the dual variable λ with corresponding subproblems. The lower boundary (LB) of (17), which is z, can be updated by finding the better upper bounds of the D(λ).
The NC in (24) is relaxed and added to the objective in (25), where H s is structured artificially for each scenario as studied in [16,29]. Therefore, dual variable λ is shared with corresponding scenarios in same nodes, assuming that the difference is indicated by H s respective to the same radius λ as depicted in Figure 2. The underlying problem becomes maximising objectives concerned with λ, andx is ignored as derived in the second line of (25). It should be noted that subproblems are still minimising objectives regarding x s as denoted in the third line of (25).
The differential of D(λ) and its expected bounds z are used to update λ. We can derive (26) by differentiating (25) with λ. The z in (27) is decision variables for the LBs for the dual function, and its differential is decided by step size p. PB uses the stability centre pointλ to represent the implementable value of λ and it decides the better LB of (17), giving higher value of bounds with the smallest possible radius from the optimal point as depicted in Figure 2. Therefore, a new λ is obtained by maximising the following optimisation problem as follows: In (28), z s is maximised to satisfy its optimal condition ∂z − p k (λ −λ) ∈ 0, which is the maximisation of the objective with constraints for its cutting planes that are indicating the possible movement for z s with the differential of subproblems (26). During iterations, k ∈ [1, ..., k max ] , PB updatesλ only if the improvement is higher than the plausible minimum value.
The differential of subproblems given in (29) represents the actual improvement of solutions with a new λ k . The difference between the feasible boundary z s and current boundary D(λ k−1 ) given in (30) represents the expected improvement of subproblem solutions with λ k . PB accepts new λ asλ when the improvement of the solution which can be confirmed by comparing subproblem solutions, ∂D(λ) to ∂z. If the improvement is higher than mL it will updateλ, or else it will go through the null step [29].
Thus, the stability centreλ can be updated if the actual improvement is larger than mL ∈ (0, 0.5] of the estimated improvement as given in (31). The convergence of variables is thus indicated when there is no expected improvement of the bounds of (27) such that ∂z ∈ 0. In each iteration, the step size p k is calculated according to results of (31) and (32) to smooth the heuristic searching process. When the actual improvement is too small, the step size decreases; on the other hand, when the actual improvement is sufficiently high with mR ∈ [0.5, 1) as given in (32), the step size increases. The step size is determined as follows: In (33), new step size for the next iteration h k is yielded by assuming that the actual improvement is the same as the expected improvement as follows: Therefore, p k+1 is adjusted by h k and initial configurations.
In (35), the step size is decreased when (31) is satisfied, and in (36), the step size is decreased when (32) is satisfied. Note that h k ≤ u k when the factor mR ∈ [0.5, 1) as follows: The detailed process to update the step size is featured in [17], and the detailed process of the PB method is summarised in Algorithm 2.

Branch and Bound in MSSP
Solutions of MSSP given in Section 3.1 can be obtained by using the methods suggested in former subsections. However, the solution is infeasible if the integrality and the NC are not satisfied. The BB algorithm bounds unconverged variables with all corresponding subproblems. A search tree is generated to observe varying results of the bounding process for the NC and manages multiple scenario trees in it; the bounding process is not limited to the single MILP problem. Moreover, the gap between the true solution of (17) and the boundary of (25) is resolved by controlling the additional cost from BB constraints. Thus, the gap between the true solution of (17) and the boundary of (25) is resolved during the process of BB.
First, after the process of Algorithm 2 is done, it generates a search tree where branches B are structured to find solutions with different conditions. Constraints added by the branching are applied to all corresponding subproblems. These constraints are added in addition to the original constraints G, so that G B includes G for all branches in the B. Second, to branch variables which violate integrality, it splits B into B and B by adding two disjunctive constraints.
In (38) and (39), x t is a node variable for time stages t = [1, ..., T], and the projection in (18) is used to generate the integrated solutionx t . We can choose a node variable having the most fractional value. On the other hand, it is possible to branch over continuous variables where the NC is not satisfied.
In (40) and (41), NC is a tolerance of the NC to make a disjunction. The measurement of violation can be used to choose the variable, e.g., the highest distance between subproblems.
Termination of BB is determined when there is no need to branch further or no feasible solution. We can define a node as fathomed when the branching on this node is clearly meaningless. This node becomes a permanent leaf. The process will be terminated when all branches are fathomed. If all variables satisfy NC and the objective value is lower than the current best upper boundary (UB) and also feasible, the node is fathomed as the best UB. On the contrary, a node can be also fathomed if the solution is infeasible or is higher than the best UB. The best LB is designated by the lowest boundary among leaf nodes that is not fathomed, which makes it possible to have better results than the current best UB. Therefore, the BB algorithm can also be terminated when LBs closely approach the best UB.
Since BB yields the higher boundary of (25), so that Z LB ends up with the same value of the UB or higher, it is obviously not required to branch on nodes that satisfy the inequality in (42). Finally, Z UB is accepted as the final solution of MSSP that obtains the condition of global optimality. The termination process can be summarised as follows:

1.
The node would be "fathomed" if new bound D(λ) is higher than the current optimal or is infeasible.

2.
The node would be updated as the Best UB if new bound D(λ) is lower than the current Best UB, and the NC is satisfied.

3.
The node would be updated as a candidate for the Best LB if new bound D(λ) is lower than the current Best UB, but the NC is still not satisfied.
The overall process of BB methods is detailed in Algorithm 3.

Simulation Results
Our simulation conditions are described, and the results are discussed in this section. The time horizon is considered up to 30 years with five years per stage; therefore, six stages are considered at maximum horizon. We use predicted data of economic growth and electrification rates to forecast the demand scenario tree, where both components are considering the electrical demand growth [30,31].

Test System Configuration
Our test bed system contains 22 load buses with 41 existing lines and 6 generators, modified from the IEEE 30-bus network. In the scenario tree, we consider 2, 4, 6 stages with 10 way, 7 way, 5 way splits for demand scenarios, respectively. The split for the scenario tree is the number of branchings from a node in the tree to the succeeding nodes.
Candidate lines for construction are marked with dotted lines in Figure 3 and featured in Table 2, so that the planner can construct lines, or else augment the capacity of existing lines. In Table 2, there are 12 candidate lines for the system network. Among candidates, there are 8 lines A-H (red) selected from the IEEE 30-bus network, and 4 lines I − L (blue) nearby bus 11, where additional load of 17 MW will be added. We assume that no more generation resources will be added, and the demand will increase every five years. That is, the number of generators is enough to cover expected demand loads, but power flows will not be delivered without expansions.

Stopping Criteria Testing
In this paper, estimations required to obtain valid results from the DD methods are reported through empirical results concerning the initial parameters of each method. Those parameters are associated with the intensity of the NC as the methods to deduce the optimal decision. We measure the method performance on the PH's step size γ, PB's NC matrix H and its step size p to verify the impact on the termination condition. Moreover, branching results from BB are reported.

Progressive Hedging
The step size γ represents the slope of the NC to update the penalty ω. A steep slope for the NC with a small γ takes less progressive steps for the gaps enclosing the optimal decisions. On the other hand, a gentle slope with a large γ more rigorously forces convergence. However, the simulation time does not directly follow the step size; in fact, subproblem solutions are more likely to cycle around the optimal solution when the step size increases. The cycling effect of decisions can cause a weak convergence rate that disturbs the integrated solution to satisfy the original constraints of the TLEP. As shown in Table 3, very small or large step sizes lead to expensive termination in our cases, requiring a lot of time to solve. Figure 4 depicts 4 cases with different γ values tested on the 4 stage 7 split TLEP problem. Cases with step size 1 and 5 successfully terminated, yet cases with 0.1 and 10 failed to reach the stopping criteria 0.0001 under the given 100 iteration limit. We observed that the simulation time to solve subproblems was longest using the largest step size and the objective function was not properly minimised. This case shows poor convergence rate with no further updates of the penalty for NC ( Table 3). The objective obtained by the smallest step size was better than others but failed to terminate. The convergence rate is not sufficient when the step size is too small or too large. Valid solutions were obtained with step size 1 and 5, however, the objective was not minimised compared to the value for step size 0.1. Hence, there was a weak convergence rate in later iterations, which is entailed by using a fixed step size and inevitable when the problem size increases. This clearly shows the requirement of parameter tunings to accelerate the process and obtain high-quality solutions, by continuously modifying the step size accordingly. However, it is difficult to designate the appropriate NC penalty for every decision variable.

Proximal Bundle
The measurements for the PB method are compared in Figure 5. Initial parameters described in Section 3.3 represent the step size p, which is proportional to the improvement of the dual function, and NC coefficient H s , which finds the gradient for cutting planes. We measured the variable convergence (nonanticipativity), simulation time, and the objective cost of the TLEP problem. Given parameters are described by per-unit values based on the objective cost and the stopping criteria was 0.0001. Given NC coefficient H s in Figure 5, different simulation results were observed according to the initial step size. When the initial step size is very low or H s is very large (grey), simulations take more steps with little improvements to the solution. Iterations wasted to draw redundant cutting planes take a long time to deduce the dual function. On the contrary, immediate terminations were observed for the highest p of 100 respective to 0.5 and 1 H s . However, they did not perform the best on the NC results since expected new boundary of the dual function was not able to be improved without changes of dual variable λ. Among results with the same H s (green, magenta), the step size 0.1 performs the best for the NC and objective. As a result, we verified that the wrong choice for parameters can bring premature termination or expensive termination. It is difficult to distinguish which cost performs the best as the process is premised to incline from the minimum to the maximum, yet the original problem involves minimisation. Thus, we excluded the high costs that violated the NC from acceptable solutions.

Branch and Bound
The BB algorithm was applied to the results from the previous subsection, and the associated results are analysed by comparing the best performance and the worst performance in the same search tree.
Intermediate n-th branching results until they are fathomed are marked in Figure 6, where the best and the worst are distinguished. The best case can quickly converge integer variables with low additional cost for the NC, which support the objective function closely approaching to the optimal value. However, the result highlights the importance of first branching nodes, and objectives increased by the NC do not always guarantee the optimal solution. In Figure 7, results for every node in the branching tree are marked, in which nodes are distinguished, following the best case or the worst case. There were significant differences between the two groups. First, the best-case nodes converge close to the optimal cost, however, the worst-case nodes dramatically increase for the first few branchings and produce an overestimated cost. Second, the worst-case nodes spend more time solving additional PB processes to converge the integer variables. Therefore, the quality of the objective function is not guaranteed if we follow the worst branching results and this is determined from almost the first branching.

Numerical Results Analysis
Based on the empirical test for the initial parameters, we report the TLEP results based on MSSP with different program sizes in Table 4. Line construction and generation cost amounts are displayed to show which line has been augmented or constructed according to DD methods. The termination conditions considered two points: method termination and feasibility. Method termination indicates whether the stopping criteria was reached within a time limit of 48 hours. Node feasibility indicates that the integrated node solutions from the projection satisfy constraints involved in the TLEP. In general cases, if obtained solutions satisfy the NC, the node feasibility also is satisfied.
The PH method contains quadratic terms in the objective function in subproblems and updates penalty parameters andx by projection. The quadratic terms in the subproblems usually take more time to solve single subproblems than that of PB. However, the PB method has the main optimisation problem of updating the dual variable λ, which takes more time than calculating the projection. According to the numerical results, solutions from the PH method satisfy the NC because its termination condition pertains to the NC. However, the PB method can violate the NC since the termination is based on boundary improvement not on variables. The DDSIP which implements the BB algorithm after the PB termination finds integer variables implementable across corresponding subproblems. Therefore, we used medium-iteration penalties and concatenated them for the warm-start of the PB method to speed up the process. We convert the penalty parameters for the PH method as follows: The parameter exchange enables the penalty matrix from the PH method to support the PB method to find the initial point. Among the cases featured in Table 4, constructions for new lines occur in Case 2 and 3, whereas Case 1 decides to upgrade existing lines. In general, construction for a new line is more expensive than upgrades so the results reveal the cost savings that construction can provide over 30 years. The selected candidate lines are J and L, which would lower the cost for generating power flow by linking end bus 11 to the central area off the network. Line 8-28 was the most frequent upgrade choice, which is the connection for end bus 8. On inspection, there is a candidate line E for bus 8 that obviously links the network to bus 8 but this is not as beneficial as linking the bus 11. As shown in Case 1, PH and PB determine an upgrade for 16-17 whereas others determine an upgrade for 21-22; the cost effectiveness can be verified by the value of objectives and proximity to the EF results. Growing the system demand, it becomes obvious that the TLEP problem in Case 2 requires new lines as the PH method gave the worst objective in spite of using the best plan from Case 1 (note that there were no more expansions after the first stage in Case 2). The expansion for J is displaced by L in Case 3, observing that most of the generators were fully contributed in a 30-year time span and subproblems determined expansion L and upgraded 21-22 to work better with other plans. This change accounts for the more severe demand conditions of the few scenarios where further upgrades to 23-24 after the first stage are decided.
Specific optimisation results are shown in columns 7-10, Table 4. The objective comprises the cost for generation, line expansions, and supply deficits. This paper depicts the sum of costs for individual results to illustrate all obtained objectives for comparison (even though the solution is infeasible). Generation cost rises as expensive fuel generators are committed, which can be seen when phase angles are not physically modifiable or line capacity is limited. Therefore, line expansions to provide energy from cheap generators can minimise the generation cost, yet the mitigated cost should be balanced against the cost of deficits.
Despite the decrease in speed as problem size is enlarged and the slight overestimation of the objective entailed by projection, PH in [14] sufficiently satisfies the NC if it can terminate and it obtains the solutions in a reasonable time. Whereas PB in [17] rarely succeeds at the same feasibility test; actually, the results of PB rank lowest but are not implementable (the expansion result in Case 1 merely fulfils the integrality). Feasibility can be satisfied by the BB algorithm with an adequately reduced objective as verified in DDSIP [16] and PH + DDSIP (proposed) results. Moreover, simulation time for the DDSIP process is lowered by the penalty matrix of the PH. The reconciliation of two different DD methods is observed from this result that mitigates both problems for finding the parameter organisation in the PB and the cycling effect that the PH undergoes.
The EF results formulated as (15) can represent a reference result for the optimality in Case 1 and 2, but it requires myriad time and a larger MIP gap. Considering the test bed is based on a simplified 30-bus network with some electrical equipment omitted, it is not suitable to use EF in a practical evaluation process. In Case 3, the EF failed to solve the problem within the given time. Instead, the objective 6256.58 was obtained through a 3-split problem, which was slightly overestimated. Considering more detailed scenario trees usually lower the cost, the expected cost for extreme scenarios is discouraged in DDSIP and PH + DDSIP for Case 3 as the gap between scenarios is smaller compared to the 3-split case. As a matter of fact, the objective does not necessarily need to be the same to obtain the unified results as shown in Case 3. DDSIP has a higher objective as many scenarios chose expansions to avoid supply deficits. On the other hand, PH + DDSIP chose expansions in only a few extreme scenarios, and the additional cost for the deficit was not that significant. Regardless of the extremes, expansions for L and upgrades for 8-28, 21-22 are necessary to provide energy throughout the 30 year time span. Particularly, as the number of scenarios increases as the scenario tree branches, so probabilities for single scenarios decrease and an influence to first stage decision is shrunk. Finally, required transmission lines for the network and times needed for completion can be estimated in chronological order.
The purpose of our simulation was not to find the one true solution for the test system. Yet, we focused on verifying variable convergence of MSSP problem, composed of typical TLEP problem, and we extended up to six stages. In this long-term evaluation process, we accomplished the feasibility of the results making the integrated variables from DD methods implementable. Various test cases revealed that the expansion plan can be relative according to method type regardless of the feasibility condition. Not only considering the accuracy, our test revealed the problem of time. The PH method shows contradictory results depending on the problem size where decision variables do not converge despite many iterations of the long-term problem. Whereas the proposed method that combines PH and DDSIP shows coherency over various sizes. Moreover, expansion results over different time spans sufficiently explain the consistency in response to the electricity demand. Many studies have focused on the financial benefit of resources in the objective function. It represents the expense, which apparently appears during the production or is induced subsequently [7,28,32]. Instead of that, we provide a framework for DD methods to handle integer variables in MSSP. Through the experiments, algorithms are tuned to cover multistaged problems and they give variety to the assessment results when there is uncertainty in the network organisation.

Conclusions
In this paper, we introduce the methodological basis of dual decomposition methods for multistage stochastic programming, tailored toward the transmission line expansion problem. Since the role of transmission lines is as substantial as generation resources are, careful planning is required as expansions are not retractable. We evaluate this investment decision on the 30-bus test network by considering three costs-reliability, operation and investment-based on future demand uncertainties. Despite the shortcoming in integer variable convergence in conventional dual decomposition algorithms, the proposed method enables the utilisation of several scenarios over a long-term horizon by finding correlation between scenario subproblems represented by the nonanticipativity constraint. Furthermore, we bring three dual decomposition methods-progressive, proximal bundle, dual decomposition in stochastic integer programming-and combined them to account for the convergence issue by adjusting warm-start values for internal parameter settings. The results not only convey optimal decisions among furcated future realisations but also achieve the economical use of time. We illustrate the significance of a long-term evaluation by showing that empirical tests over different time spans can yield various combinations of expansions, enabling planners to confirm them according to potential benefits under uncertainty. However, further studies for realising adequate reliability cost and evaluating power system flexibility with newly constructed lines should be considered in the future studies.