A MEAN-VARIANCE SERIAL REPLACEMENT DECISION MODEL: THE CORRELATED CASE Matt Brown Department of Industrial Engineering and Operations Research Etcheverry Hall University of California, Berkeley Berkeley, CA 94720 Technical Report 92-10 January 1992

A Mean-Variance Serial Replacement Decision Model: The Correlated Case Matt Brown Department of Industrial Engineering and Operations Research Etcheverry Hall University of California, Berkeley Berkeley, CA 94720 (510) 643-7739 Internet: mjb@um.cc.umich.edu Technical Report, Department of Industrial and Operations Engineering, College of Engineering, The University of Michigan, Ann Arbor, MI (January 1992).

A Mean-Variance Serial Replacement Decision Model: The Correlated Case Abstract An optimal approach for solving serial replacement problems with uncertain rewards and a risk-averse decision maker was previously introduced by Brown. In this note, we extend Brown's model to allow rewards to be correlated across replacements. Correlation may be useful in modeling the effects of the external environment on the rewards. However, total enumeration appears to be the only method to find optimal solutions when rewards are correlated. Thus, we develop a heuristic and compare it to several existing solution approaches. Computational results indicate that, on the average, our heuristic outperforms the other solution approaches by at least 8 percent for positive correlation and 17 percent for negative correlation. KEY WORDS: MACHINE INVESTMENT; STOCHASTIC NETWORKS: LONGEST PATH PROBLEMS; DECISION MAKING UNDER UNCERTAINTY. i

A Mean-Variance Serial Replacement Decision Model: The Correlated Case 1 Introduction A dynamic programming procedure to find all Pareto-optimal solutions for serial replacement problems with uncertain rewards and a risk-averse decision maker was previously introduced by Brown [3]. A solution specifies the sequence of assets that are replaced over time, including the installation and replacement times of each asset. The service life of an asset is the time that elapses between installation and replacement. In this note, we extend Brown's model to allow rewards to be correlated across replacements, referring to this as the correlated case (as opposed to independent case). We emphasize that the correlation being modeled is across replacements, and we differentiate this type of correlation from that across different service lives for the same asset (e.g., the reward for an asset providing three periods of service being correlated with the reward for the same asset providing four periods of service). Correlation of this latter type can be included in both Brown's model and our extension. Correlation across replacements may be useful in modeling the effects of the external environment on the rewards. For example, consider a decision maker forecasting the operating expenses of a delivery vehicle, who is uncertain about the route and associated road surface conditions for the vehicle. To reflect this uncertainty, the decision maker specifies a range to cover all reasonable outcomes. However, once the vehicle is placed in service, the decision maker expects the route and conditions to 1

remain relatively constant over time. Thus, the operating expenses of a future vehicle will be correlated with the current vehicle that it will replace. Our investigation has two research objectives. First, we want to develop a solution procedure. When correlation is considered, the resulting problem structure is such that total enumeration appears to be the only optimal solution approach. Thus, our solution procedure is a heuristic. Second, we want to compare our heuristic both to an upper bound and to several existing solution approaches. The remainder of this paper is organized as follows. In Section 2, we state the assumptions of our model. We describe our heuristic and two different error bounding procedures in Section 3. In Section 4, we discuss our computational study, and summarize our findings in Section 5. 2 Assumptions of the MV Decision Model The major assumptions of the MV decision model are as follows. * Assets are replaced serially. * The decision maker is uncertain about the monetary reward that would result from installing and operating an asset. * The reward for the service provided by an asset is summarized by its after-tax net present value (NPV) distribution. These rewards can be either positive or negative. 2

* The NPV of an asset being appended to a partial sequence is correlated only with the NPV of the immediately preceding asset in the sequence. Correlation is limited to first order effects because it captures the major effects and a more elaborate model would require a decision maker to forecast all the higher orders of correlation, which may quite difficult. * The NPV of an asset is normally distributed, and hence the NPV of any sequence of assets will also be normally distributed. * Each asset has a known, maximum service life beyond which it must be replaced. * There is a finite number of alternative assets from which to choose at any point in time. * The decision maker's utility function is a continuous, concave, monotonically increasing function of money (w). * The decision maker's objective is to maximize expected utility. The problem then is find the sequence of assets of highest expected utility (EU = f U(w)f(w)dw, where U(w) is the decision maker's utility function and f(w) is the normal density function of the NPV for a sequence of assets) over all possible sequences that provide service from time zero to the finite horizon (H). 3 Solution Approaches When rewards are correlated across replacements the dynamic programming principal of optimality does not hold because the optimal solution does not necessarily consist 3

of optimal partial solutions. We provide an example below. Consequently, we cannot identify optimal solutions using dynamic programming. The only method to find optimal solutions appears to be total enumeration (see Bard and Bennett [1]; and Loui [4]), and the associated computational burden can be prohibitive. Consider a problem with five asset types available at each point in time, and with each asset having a maximum service life of just one period. For a planning horizon of 20 periods, the total number of solutions is 9.54E + 13. If the maximum service life increases to two periods, the total number of solutions increases to 1.95E + 15. ( Given the difficultly of obtaining optimal solutions for the correlated case, we develop a heuristic. Our heuristic is an extension of the dynamic programming approach Brown [3] uses for the independent case. In the remainder of this section, we summarize Brown's approach, then we present our heuristic, and finally we describe the two error bounding procedures that we use in our computational study. 3.1 Brown's Approach for the Independent Case The approach of Brown [3] uses dynamic programming to find all Pareto-optimal solutions (sequences of assets). At each time period (stage), a Pareto-optimal set of sequences is identified. The Pareto-optimal set of sequences also is a MV-efficient set of sequences, under the assumption that the NPV of each asset in a sequence is normally distributed. Letting St be the set of all sequences providing service from time zero to time t, a sequence s E St with mean Mu(s) and variance u2(s) is MV 4

efficient if for every s' E St, s' = s, u(S) > p(S') and a2(s) = a2(s'), or,(s) >,u(s') and a2(s) < 2(s'). The dynamic programming algorithm proceeds in a forward direction, with the stages being the time periods and the states being the asset type (j) and time of installation (t). The objective function values are the NPV mean and variance of a given sequence. At each point in time, the decisions are what asset type to install and for what length of service. Let S7 C St be the set of all MV-efficient sequences providing service from time zero to time t, Jt be the set of assets available at time t, and Nt(j) be the maximum service life for asset j available at time t. The set St can be formed as follows St = {sU (j,t' - t'): s E St,j E Jt, t' > 0 t - t' = 1,2,..., Nt,j) where (j, t', t - t') represents asset type j, installed at time t', for t - t' periods of service. Expressed differently, St is the set of all sequences formed by appending one additional asset, j E Jtt, to the end of a MV-efficient sequence providing t', t' < t, periods of service, such that the resulting sequence provides t periods of service. St is then given by S; = {s E St: s is MV-efficient}. The boundary condition is So = 0. Brown [3] notes that some problems cannot not be solved optimally due the state space becoming prohibitively large. In these cases, Brown uses a heuristic called the 5

cluster heuristic. The cluster heuristic identifies a subset of Pareto-optimal solutions by keeping only those sequences that are not stochastically dominated if the NPV distributions are truncated 6 standard deviations from the mean, where 6 is a parameter adjusted by the heuristic. 3.2 A Heuristic for the Correlated Case The heuristic we develop for the correlated case is a slight modification of Brown's optimal approach. The modification involves adding the covariance term to the variance. For example, if j is the asset to be appended to a partial sequence and it follows asset i, then the variance term is given by j + 2pijaiaj, where pij is the correlation coefficient for j following i. The approach is a heuristic because it will fail to identify the optimal solution if the optimal path is not MV-efficient at each replacement epoch. To illustrate why, consider the following example with a horizon of two periods. At time zero, two assets, i and j, are available, both with identical one period service lives. At time one, only asset k is available and only for one period of service. Let pi =j, a2 =64, j = 81, a,2 = 100, Pik = 0.5, and Pjk = 0.1. For this example, there are two possible sequences: (i) asset i followed by asset k, [i, k], and (ii) asset j followed by asset k, [j, k]. In comparing these two sequences the means are identical and thus can be ignored. For [i, k], the sequence variance is 64 + 100 + 2(0.5)(8)(10) = 244. For [j, k], the sequence variance is 81 + 100 + 2(0.1)(9)(10) = 199. Thus, the optimal sequence is [j,k]. However, our approach will select [i, k] because at time 1, ([i, a2) MV-dominates (I,j,?a). UL~trU L, VI11~,$ E- Z )V jj 6

We also modify Brown's cluster heuristic to use for problems in which the state space becomes prohibitively large. Again, the modification involves adding the covariance term to the variance. Otherwise, the procedure is exactly the same, including the method for selecting which sequences to eliminate. We report the number of times we use the cluster heuristic modification for the problem sets considered in our computational study. 3.3 Bounding Procedures Given total enumeration appears to be the only method to identify optimal solutions and such an approach can be computationally prohibitive, we decided to compare the solutions obtained by our heuristic to an upper bound. However, obtaining a tight bound on the optimal solution is also a difficult problem. Consequently, we use two different bounding procedures in our computational study, selecting the smaller (tighter) of the two bounds. The first procedure computes a bound using the maximum mean and a lower bound on the variance across all sequences. We will refer to this procedure as the minimum variance bounding procedure. The maximum mean and minimum variance will MV-dominate all other sequences and hence the resulting expected utility of the sequence will be an upper bound. To find the maximum mean, we solve a longest path problem using dynamic programming to find the sequence of assets with the highest mean NPV. We find a lower bound on the variance by solving a shortest path problem, again using dynamic programming. The dynamic program proceeds in a forward direction, 7

with the stages being the time periods, and the states being the last asset type in the partial sequence and the last asset's minimum and maximum standard deviation for any service life. The optimal value function, SeqVar*(t), is a lower bound on the variance for any sequence providing service for the first t periods. The key is to find a tight lower bound on the covariance. To do this, at each stage for each asset type j, a lower bound on the variance is found for any sequence in which the last asset is of type j. Let SeqVar*(j, t) be this lower bound, StdMin(j, t) be the minimum possible and StdMax(j, t) be the maximum possible standard deviation for asset type j in any sequence providing t periods of service for which j is the last asset, and o(j, t', t - t') be the standard deviation for asset type j, installed at time t', providing t - t' periods periods of service. Then: StdMiin(j, t) = Min {a(j, t', t - t') t - t'= 1,2,..., Nt,(j), t' > 0} StdMax(j, t) = Max{((j, t', t - t'): t ' - 1, 2,..., N,(j), t' > 0} and the functional equation is SeqVar*(j, t) = Min {SeqVar*(j', t') + a2(j, t', t - t') + 2pj,jStdMin(j', t')a(j, t', t -t'), SeqVar*(j', t') + a2(j, t', t - t') + 2pjjStdMax(j', t')a(j, t', t - t'): t -t' = 1,2,..., Nt(j), t' > O, j E Jt} At each stage, SeqVar*(t) = Minj{SeqVar*(j, t)}. At the horizon, if SeqVar*(H) < 0, a minimum variance of zero is used. The boundary conditions are SeqVar*(j, 0) 8

= for all j. The second bounding procedure is only used for problems in which all the correlation coefficients are nonnegative. In this case, the variance of any given sequence assuming independence will always be less than or equal to the actual variance, and the associated expected utility for the independent sequence will always be greater than or equal to the expected utility for the same sequence when the covariance is included. Thus, an upper bound can be found by solving the problem assuming independence. We will refer to this procedure as the independent bounding procedure. 4 Computational Experiments The goals of our computational experiments are (i) to evaluate the error bound performance of our heuristic, and (ii) compare the sequence identified by our heuristic to the sequences identified using three existing solution approaches. (The same three Brown [3] uses for the independent case.) These approaches are the traditional (TRAD) approach, the expected value (EV) approach, and certain monetary equivalent (CME) approach. Each is designed to solve a problem that is related to, but not the same, as our problem. We distinguish among the different sequences by using superscripts (e.gv., sTRAD for the sequence identified by the traditional approach). We denote the sequence identified by our heuristic as sUTIL.e The traditional (TRAD) approach assumes an infinite horizon and that all future assets are monetarily identical to those available currently (see Newnan [5]). This approach replaces each random NPV by its mean and then selects the asset type and 9

service life pair with the highest expected annual equivalent value (NPV amortized over the service life). This decision is then repeated forever. The expected value (EV) approach substitutes the mean NPV in place of the NPV distribution, and then finds the sequence of assets with the highest mean NPV. Oakford, Lohmann, and Salazar [6] present a dynamic programming algorithm that can be used to find such a sequence. The certain monetary equivalent (CME) approach replaces each random NPV by its certain monetary equivalent. The certain monetary equivalent is the constant (k) yielding the same utility as the random NPV (i.e., U(k) = EU). The solution is then obtained by dynamic programming in exactly the same manner as for the expected value approach. In the remainder of this section, we first define the performance measures for our study. Then, we discuss how our problem sets are generated. We note that the performance measures and test problem generation process are identical to those used by Brown [3] for the independent case, except that we add a procedure to generate the correlation coefficients. Finally we report computational results which indicate that, on the average, our heuristic outperforms the other three solution approaches by at least 8 percent for positive correlation and 17 percent for negative correlation. 4.1 Performance Measures We collected four performance measures: (i) the CPU time required by each approach, (ii) the maximum cardinality of any MV-efficient set of sequences at any stage required 10

by the optimal approach or cluster heuristic, (iii) the percentage of sequences found by a given approach that match (have expected utility equal to) sUTIL, which we will refer to as the matching performance, and (iv) the utility performance of a sequence compared to sUTIL, which is given by EU(s) - EU(sRAND) EU(sUTIL) -EU(SRAND) where sRAND denotes the sequence of highest expected utility out of 100 random sequences. The fraction is set to one if EU(sUTIL) = EU(sRAND), and to zero if EU(s) < EU(sRAND). We note that the most obvious utility performance measure, EU(s) / EU(sUTIL), cannot be used. At best, this fraction is difficult to interpret, and it may even be undefined because both the NPV and expected utility of a sequence may be negative, zero, or positive. By using a random sequence as a benchmark of performance, the resulting performance measure ranges from zero to one. Random sequences are formed by randomly selecting the asset type and service life at each replacement epoch. For our heuristic, the last two measures are computed by comparing sUTIL to an upper bound. The resulting utility performance is given by EU(sUTIL) _ EU(sRAND) EU(sUB) - EU(sRAND) where sUB is the sequence of highest expected utility found by the upper bounding procedure. The value is set to zero if EU(sUTIL) < EU(sRAND). 11

Table 1: Parameter Ranges Parameter Low High Discount Rate U(0.1,0.2) U(0.2,0.3) Horizon DU(10,25) DU(25,40) Service life DU(2,8) DU(8,14) Coefficient of Variation U(0.1,0.4) U(0.6,1.2) Difference U(0.0,0.05) U(0.05,0.1) Risk Aversion U(1.5,5) U(16.5,20) U - uniform distribution, DU _ discrete uniform distribution 4.2 Test Problem Generation To generate a broad range of problems, we use a 26 factorial design. The six factors and their ranges are listed in Table 1. We first draw values for the discount rate, the horizon, the number of asset types available at any point in time (which ranges from two to seven for our problem sets), and the maximum service life for each asset type. Second, we generate the NPV means and variances for every possible service life for each asset type, using the coefficient of variation and difference factors. Third, we set the decision maker's risk aversion, continuing to use the same three mathematical forms for the utility function that are used in the independent case: (i) the exponential, U(w) = (1 - ecw)/c, reflecting constant risk aversion, r(w) = c, (ii) the logarithmic, U(w) = ln(w + b), (w + b) > 0, capturing decreasing risk aversion, r(w) = 1/(w + b), and (iii) the power, U(w) = (w - wo)3, wo < w for all w, 0 < p < 1, also reflecting decreasing risk aversion, r(w) = (1 - P)/(w - wo), where r(w) is the risk aversion function. Finally, we generate the correlation coefficients. We restrict the correlation coefficients to be equal across all possible service lives for a given asset type (i.e., the correlation coefficient for asset j following asset i, 12

pij, is independent of the service life of either i or j) to simplify the generation process. However, note that in practice our heuristic can be used for any level of detail specified by the decision maker. In modeling correlation across replacements, positive correlation typically arises more frequently than negative correlation (since situations in which there is some substitution effect are less common). Consequently, we study two cases: positive correlation and negative correlation. For positive correlation, all the coefficients generated are nonnegative. Let u(a, b) represent a random draw from a uniform distribution with endpoints a and b. The coefficient for sequential assets of the same type is set equal to u(O, 1). The remaining values for the coefficients for the current asset type followed by an asset of a different type are set by multiplying u(0, 1) times the coefficient just set for sequential assets of the same type. For negative correlation, the coefficient for sequential assets of the same type is set equal to zero. The remaining values for the coefficients for the current asset type followed by an asset of a different type are set equal to u(O, 0.5), with the sign set based on the observation that all the coefficients cannot logically be negative (for more than two assets). As an illustration, consider three assets. If the first and second assets are negatively correlated, and the second and third assets negatively correlated, then it does not make sense for the first and third assets to also be negatively correlated. Therefore, the negative case is actually a mixture of positive and negative values. 13

Table 2: Utility Performance by Utility Function-Positive Correlation Utility Decision Percent Utility Performance Function Procedure Matching Average Minimum Maximum Std. Dev. Exponential TRAD 34.37% 0.5308 0.00 1.00 0.46 EV 63.75% 0.7720 0.00 1.00 0.39 CME 84.06% 0.9230 0.00 1.00 0.26 UTIL 1.56% 0.7428 0.00 1.00 0.29 Logarithmic TRAD 35.00% 0.5557 0.00 1.00 0.46 EV 66.25% 0.8079 0.00 1.00 0.37 CME 39.06% 0.7079 0.00 1.00 0.40 UTIL 3.43% 0.7740 0.00 1.00 0.28 Power TRAD 34.68% 0.5593 0.00 1.00 0.46 EV 65.00% 0.8058 0.00 1.00 0.37 CME 31.87% 0.6073 0.00 1.00 0.46 UTIL 3.43% 0.7717 0.00 1.00 0.28 4.3 Results We randomly generated 640 problems-320 with positive correlation and 320 with negative correlation-by replicating our 26 design five times. For positive correlation, 118 problems required the cluster heuristic modification, while for negative correlation, 141 problems required it. Broadly speaking, these results indicate problems with negative correlation are more difficult to solve, primarily because negative correlation drives the variance towards zero which makes it more difficult to distinguish among the different assets. Both are higher than the 92 times the cluster heuristic was required for the independent case. Summary statistics for the matching and utility performance measures appear in Tables 2 and 3. The third column shows the average matching performance. The last four columns show the utility performance average, minimum, maximum, and standard deviation. While it may appear that our heuristic (UTIL) is outperformed by at least one other approach for each utility 14

Table 3: Utility Performance by Utility Function-Negative Correlation Utility Decision Percent Utility Performance Function Procedure Matching Average Minimum Maximum Std. Dev. Exponential TRAD 35.00% 0.5460 0.00 1.00 0.46 EV 66.56% 0.7966 0.00 1.00 0.37 CME 79.06% 0.8942 0.00 1.00 0.29 UTIL 4.06% 0.6051 0.00 1.00 0.35 Logarithmic TRAD 35.93% 0.5670 0.00 1.00 0.46 EV 70.00% 0.8381 0.00 1.00 0.34 CME 27.81% 0.6682 0.00 1.00 0.39 UTIL 4.06% 0.6193 0.00 1.00 0.34 Power TRAD 36.25% 0.5666 0.00 1.00 0.46 EV 69.37% 0.8367 0.00 1.00 0.35 CME 21.56% 0.5827 0.00 1.00 0.45 UTIL 4.06% 0.6201 0.00 1.00 0.36 function, in fact our heuristic outperforms them all. Recall that both the matching and utility performance use the UTIL sequence as the basis of comparison when computing statistics for the TRAD, EV, and CME approaches. In this way the measures will tend to highlight any differences in the sequences found. If however, an upper bound is used as the basis for comparison (as is the case when the UTIL performance measures are calculated), the difference across the solution approaches would be more difficult to recognize unless the bound is tight. The average utility performance for the UTIL sequence is lower for negative correlation than for positive correlation. Generating tight bounds is especially difficult in the former case. However, since positive correlation will typically arise more frequently, the importance of obtaining tight bounds is reduced. Tighter bounds are obtained for positive correlation because the independent bounding procedure can also be used. For our problem set with positive correlation, the independent bound 15

ing procedure produces a tighter bound 77.2 percent of the time for the exponential utility function, 76.9 percent of the time for the logarithmic function, and 68.8 percent of the time for the power function. However for negative correlation, only the minimum variance bounding procedure can be used, and the resulting variance can be zero. In fact, for our problem set with negative correlation, the upper bound was computed using a variance of zero (i.e., at the horizon SeqVar*(H) < 0) 18.8 percent of the time for the exponential utility function, and 15.9 percent of time for the both the logarithmic and power functions. Because of the negative coefficients, the lower bound for the minimum variance is driven down towards zero, potentially underestimating the true variance by a significant amount and resulting in a looser bound. The TRAD and EV approaches require 0.3 CPU seconds, on the average (using an Apollo DSP-4000 workstation), to identify a solution. The CME approach requires about 3 CPU seconds, while the random approach requires about 8 CPU seconds (due to the fact that 100 sequences are generated and the expected utility for each computed using numerical integration). The average CPU time for our heuristic is about 9 CPU seconds, and no problem requires more than 150 CPU seconds. These results suggest that our heuristic can solve a wide range of problems in a computationally efficient manner. An analysis of variance of the CPU times required by our heuristic shows that, for both positive and negative correlation, two factors are significant (at a = 0.10): (i) the horizon, and (ii) the maximum asset service life. More CPU time is required 16

when the horizon is high or the maximum asset service life is high. These findings are identical to the independent case. 5 Conclusion The computational results indicate our heuristic can solve a wide range of correlated problems in a very reasonable amount of CPU time. The results also indicate that, on the average, our heuristic outperforms the other solution approaches by at least 8 percent for positive correlation and 17 percent for negative correlation. However, the average utility performance of our heuristic as compared to an upper bound is not especially good, and is worse for negative correlation. We believe this is a result of the upper bounds being loose, and thus consider our results to be very conservative. An interesting question is whether a bounding procedure can be developed that will yield tighter bounds. References [1] Bard, J.F., and J.E. Bennett, "Arc Reduction and Path Preference in Stochastic Acyclic Networks," Management Science, Vol. 37, No. 2, 198-215 (February 1991). [2] Brown, M.J., A Mean-Vaiance Serial Replacement Decision Model, PhD Thesis, Department of Industrial and Operations Engineering, College of Engineering, The University of Michigan, Ann Arbor, MI, 1991. 17

[3] Brown, M.J., "A Mean-Variance Serial Replacement Decision Model: The Independent Case," Technical Report #91-38, Department of Industrial and Operations Engineering, College of Engineering, The University of Michigan, Ann Arbor, MI (December 1991). [4] Loui, R.P., "Optimal Paths in Graphs with Stochastic or Multidimensional Weights," Communications of the ACM, Vol. 26, No. 9, 670-676 (September 1983). [5] Newnan, D.G., Engineering Economic Analysis, Engineering Press, Fourth Edition, 1991. [6] Oakford, R.V., Lohmann, J.R., and Salazar, A., "A Dynamic Replacement Economy Decision Model," HIE Transactions, Vol. 16, No. 1, 65-72 (March 1984). 18