Internal Rate of Return

Transcribed from the printed issue; not yet reviewed. Read it in the PDF of the issue, page 26

I hope you are all still enjoying these Zark reprints! Once again, we include a new crossword for your enjoyment and another classic APL problem. It would be interesting to see some new solutions to these problems, perhaps using Dyadic’s dynamic functions, or some J or even K. Solutions will be printed in the next issue.

1. Utility Corner: Internal Rate of Return¶

(The purpose of this column is to make you more productive by introducing you to utility functions. Think of utility functions as APL functions that have names instead of symbols. By expanding your function vocabulary, you’ll be able to write APL code that’s more concise, more efficient, and more readable.)

The term "Internal Rate of Return" (or IROR) is used to refer to the implied interest rate in effect from a given set of cash flows. Consider, for example, the following set of annual cash flows:

¯100 6 6 6 6 106

At the start of the first year, 100 is paid out. At the end of the next four years, 6 is paid in. 106 is paid in at the end of the fifth year. The IROR of this set of flows is 6% per year. To see this intuitively, view the cash outflow (¯100) as a deposit into a bank paying 6% per year. View the next four inflows of 6 as interest withdrawals at the end of subsequent years. View the final inflow of 106 as an interest withdrawal of 6 plus the return of the original 100 deposit.

To quantitatively verify that 6% is the IROR, the cash flows should sum to 0 when discounted to their present value or when accumulated to their future value using this rate:

      +/¯100 6 6 6 6 106×1.06*-⍳6
¯2.842170943E¯14
      +/¯100 6 6 6 6 106×1.06*⌽⍳6
¯4.263256415E¯14

Close enough to zero. Your task is to write a monadic function IROR that takes a vector of cash flows as its argument, one number per time period, and returns the internal rate of return, like so,

      IROR ¯100 6 6 6 6 106
.06

Our aim was to find the fastest possible IROR algorithm. The results surprised us. The most common method wasn’t the fastest.

From a mathematical perspective, the goal is to solve for the variable i such that the value of the following formula is zero for a given vector of cash flows CF. (The symbols "f(i)" are read, "a function in terms of i.")

f(i):     +/CF×(1+i)*-⍳⍴CF

Since this problem defies a direct algebraic solution, the only alternative is trial and error (called, "the method of successive approximations"). Try a value of i, say .08, in the expression and see how close the result is to zero. Then try another value, say 0.9, and see whether its result is closer to zero of further away. As you try more and more values, you’ll get better at guessing values of i that return a value of f(i) that is close to zero.

In fact, you may find it helpful to plot the results. Plot the trial values of i along the X-axis, and the result from the expression along the Y-axis:

f(i)
  |
  |            ■
  |      ■
--+-------------- i
  |  ■
  |

Graphically, f(i) describes a curve. The trial values of i and f(i) are points along that curve. From the graphic perspective, your job is to find where the curve intersects the i-axis, that is, to find the value of i where f(i)=0.

Four different approaches were suggested:

1. Extrapolate/interpolate, using a straight line f(i)=a+(b×i) fitted to the two best points so far. The "best" points are determined to be those closest to the i-axis, i.e. for which their values of f(i) are nearest zero. Draw a straight line through those points and determine the value of i where the line crosses the i-axis. That value of i becomes the next guess. Naturally, you need to make two "random" guesses to start the process.
2. Extrapolate/interpolate, using a quadratic curve f(i)=a+(b×i)+(c×i*2) fitted to the three best points. Just as the fitting of a straight line requires tow points, the fitting of a second degree curve (parabola) requires three points. The primitive APL function ⌹ simplifies the task of fitting the curve. Since the parabola crosses the i-axis at two points, only the solution nearest the other guesses is considered.

3. Extrapolate/interpolate, using a straight line that passes through the best point and has the same slope as the f(i) curve at that point. The slope of the f(i) curve at any point is the value of the derivative of f(i) with respect to i at that point. Computing the derivative of f(i) is a straightforward calculus operation:

   f'(i):    +/CF×(-⍳⍴CF)+(1+i)*-1+⍳⍴CF
   

   Given i as the best guess so far, and f(i) as the function’s value at that value of i, a small amount of algebra produces a formula for the value of i where the straight line crosses the i-axis:

   i':       i-(f(i)÷f'(i))
   

   This value of i is the next guess. This method of successive approximations is known as the Newton-Raphson method, and is the most common technique employed for solving the IROR problem. It requires only one guess to get started. 4. Extrapolate/interpolate, using a quadratic curve that passes through the best point and has the same slope as the f(i) curve at that point, and has the same rate of change in slope (second derivative) as the f(i) curve at that point. This approach extends the Newton-Raphson method to consider the second derivative as well as the first:

   f''(i):   +/CF×(-⍳⍴CF)×(-1+⍳⍴CF)×(1+i)*-2+⍳⍴CF
   

   The subsequent algebra is considerably more complex, but the resulting curve does a better job of extrapolating or interpolating. Again, only one guess is needed to get started.

Each of these four approaches successively determines values of i that are better than the ones before. You stop the iterative process once you find a value of i for which the value of the function f(i) is "adequately" close to zero, say plus or minus 1E¯11.

Which of these approaches is the fastest?

The advantage of approaches 3 and 4 is that they require a single guess to get started. The disadvantage is that the derivative (and second derivative) must be computed for each trial of i. The computation of the derivative ( and second derivative) is more involved than the computation of the function. Fewer trials are needed, but each trial is slower.

The advantage of approaches 2 and 4 is that they use a parabola rather than a straight line when interpolating. Since a straight line is a special case of a parabola, the parabola can do no worse than the straight line and, in practice, requires about half as many trials. The disadvantage is that the calculations required at each iteration are more complex, and hence slower.

To find out which approach was fastest, we timed them on a 360-element vector of cash flows. In the process of timing, we noticed that the majority of the processing time was consumed by exponentiation — even with a math (floating point) coprocessor installed. To eliminate the exponentiation, we took advantage of the fact that the exponentiation was being done on successive integral powers, and replaced the exponentiation by cumulative products. To see how this is possible, note that the following two expressions produce equivalent results:

      (÷1+i)*⍳360    ⍝ (origin 1)

      ×\360⍴÷1+i

Eliminating the exponentiation from all four approaches improved their speeds considerably. Here are the results of the timings (scaling the results so the fastest algorithm ran in 100 "units" of time):

You can see the simplest algorithm was the fastest! Also, the most common algorithm, IROR3 using * (Newton-Raphson), was the slowest. So much for mathematical sophistication.

As always, timings are affected by your hardware and your implementation of APL. You may want to time these functions on your own machine.

The functions used in the timings are listed below. In each function, the logic for both the exponential and cumulative product algorithms is included, though the exponential logic is commented out.

    ∇ R←IROR1 V;G;J;N;P;T;⎕IO
[1]    ⍝ Given a vector of cash flows V, one
[2]    ⍝ number per time period, returns the
[3]    ⍝ internal rate of return. Uses
[4]    ⍝ straight line interpolation.
[5]     ⎕IO←0
[6]     G←1.05 1.1 ⍝ Init guesses (5%,10%)
[7]    ⍝ T←-⍳⍴V ⍝ Periods for present values
[8]    ⍝ P←(G∘.*T)+.×V ⍝ Present values
[9]     N←⍴V ⍝ No. periods for present vals
[10]    P←V+.××⍀(N,2)⍴÷G ⍝ Present values
[11]   L1:R←(-/P×⌽G)÷-/P ⍝ Interp for rate
[12]    →(∨/1E¯11>|R-G)⍴L2 ⍝ Exit if ~ same
[13]   ⍝ P←P,V+.×R*T ⍝ Present val for rate
[14]    P←P,V+.××\N⍴÷R ⍝ Pres. val for rate
[15]    P←P[J←2↑⍋|P] ⍝ Pick 2 nearest zero
[16]    G←(G,R)[J] ⍝ ... & associated rates
[17]    →L1 ⍝ Try again
[18]   L2:R←R-1 ⍝ Convert to interest rate
    ∇

    ∇ R←IROR2 V;A;B;C;G;J;N;P;T;⎕IO
[1]    ⍝ Given a vector of cash flows V, one
[2]    ⍝ number per time period, returns the
[3]    ⍝ internal rate of return. Uses
[4]    ⍝ quadratic interpolation.
[5]     ⎕IO←0
[6]     G←1.075 1.05 1.1 ⍝ Guess 7.5,5,10%
[7]    ⍝ T←-⍳⍴V ⍝ Periods for present values
[8]    ⍝ P←(G∘.*T)+.×V ⍝ Present values
[9]     N←⍴V ⍝ No. periods for present vals
[10]    P←V+.××⍀(N,3)⍴÷G ⍝ Present values
[11]   ⍝ Coefficients of quadratic:
[12]   L1:C←P⌹(G×G),G,[0.5]1
[13]   ⍝ NOTE ⌹ sometimes fails on DOM ERROR
[14]    A←C[0]
[15]    B←C[1]
[16]    C←C[2]
[17]   ⍝ Get roots from quadratic formula:
[18]    R←((-B)+1 ¯1×(0⌈(B×B)-4×A×C)*0.5)÷A+A
[19]    R←R[>/|R-G[0]] ⍝ Pick the best root
[20]    →(∨/1E¯11>|R-G)⍴L2 ⍝ Exit if ~ same
[21]   ⍝ P←P,V+.×R*T ⍝ Present val for rate
[22]    P←P,V+.××\N⍴÷R ⍝ Pres. val for rate
[23]    P←P[J←3↑⍋|P] ⍝ Pick 3 nearest zero
[24]    G←(G,R)[J] ⍝ ... & associated rates
[25]    →L1 ⍝ Try again
[26]   L2:R←R-1 ⍝ Convert to interest rate
    ∇

    ∇ R←IROR3 V;D;F;G;N;T;T1;U;V1;⎕IO
[1]    ⍝ Given a vector of cash flows V, one
[2]    ⍝ number per time period, returns the
[3]    ⍝ internal rate of return. Newton-
[4]    ⍝ Raphson: for a function f(x),
[5]    ⍝ whose derivative is f'(x), the next
[6]    ⍝ guess should be x-f(x)÷f'(x). That
[7]    ⍝ is, use the slope at x to extrap.
[8]    ⍝ to where f(x) is zero.
[9]    ⍝     f(x) = +/(G*T)×V
[10]   ⍝    f'(x) = +/T×(G*T-1)×V
[11]    ⎕IO←0
[12]   ⍝ T←-⍳⍴V ⍝ Periods for present values
[13]   ⍝ V1←T×V ⍝ Product for derivitave
[14]   ⍝ T←T-1 ⍝ Exponents for derivative
[15]    N←⍴V ⍝ No. periods for present vals
[16]    T1←¯1-⍳N ⍝ Product for derivative
[17]    G←1.1 ⍝ Initial guess (10%)
[18]   L1: ⍝ D←V1+.×U←G*T ⍝ f'(x)
[19]   ⍝ F←V+.ׯ1↓1,U ⍝ f(x)
[20]    F←+/U←V××\N⍴÷G ⍝ f(x)
[21]    D←T1+.×U÷G ⍝ f'(x)
[22]    R←G-U←F÷D ⍝ New rate
[23]    →(1E¯11>|U)⍴L2 ⍝ Exit if ~ same
[24]    G←R ⍝ Make this the next guess
[25]    →L1 ⍝ Try again
[26]   L2:R←R-1 ⍝ Convert to interest rate
    ∇

    ∇ R←IROR4 V;A;B;C;D;F;G;N;S;T;T1;T2;U;V1;V2;⎕IO
[1]    ⍝ Given a vector of cash flows V, one
[2]    ⍝ number per time period, returns the
[3]    ⍝ internal rate of return. Newton-
[4]    ⍝ Raphson extended to the 2nd deriv.
[5]    ⍝ For a function f(x), whose deriv.
[6]    ⍝ is f'(x) and 2nd deriv. is f"(x),
[7]    ⍝ the next guess should be the root
[8]    ⍝ of the quadratic equation
[9]    ⍝     0 = f(Z) = (A×Z*2) + (B×Z) + C
[10]   ⍝ where:
[11]   ⍝     A = f"(x)÷2
[12]   ⍝     B = f'(x)-f"(x)×x
[13]   ⍝     C = f(x)+(-x×f'(x))+(x*2)×f"(x)
[14]   ⍝ and:
[15]   ⍝      f(x) = +/(G*T)×V
[16]   ⍝     f'(x) = +/T×(G*T-1)×V
[17]   ⍝     f"(x) = +/T×(T-1)×(G*T-2)×V
[18]    ⎕IO←0
[19]   ⍝ T←-⍳⍴V ⍝ Periods for present values
[20]   ⍝ V1←T×V ⍝ Product for derivitave
[21]   ⍝ T←T-1 ⍝ Exponents for derivative
[22]   ⍝ V2←T×V1 ⍝ Product for 2nd deriv
[23]   ⍝ T←T-1 ⍝ Exponents for 2nd deriv
[24]    N←⍴V ⍝ No. periods for present vals
[25]    T1←¯1-⍳N ⍝ Product for derivative
[26]    T2←T1-1 ⍝ Product for 2nd deriv
[27]    G←1.1 ⍝ Initial guess (10%)
[28]   L1: ⍝ S←V2+.×U←G*T ⍝ f"(x)
[29]   ⍝ D←V1+.×U←¯1↓(÷G),U ⍝ f'(x)
[30]   ⍝ F←V+.ׯ1↓1,U ⍝ f(x)
[31]    F←+/U←V××\N⍴÷G ⍝ f(x)
[32]    D←+/U←T1×U÷G ⍝ f'(x)
[33]    S←T2+.×U÷G ⍝ f"(x)
[34]    A←S÷2
[35]    B←D-S×G
[36]    C←F-G×D-G×A
[37]   ⍝ Get roots from quadratic formula:
[38]    R←((-B)+1 ¯1×(0⌈(B×B)-4×A×C)*0.5)÷A+A
[39]    R←R[>/|R-G] ⍝ Pick best root
[40]    →(0.00000000001>|R-G)⍴L2 ⍝ Exit if ~ same
[41]    G←R ⍝ Make this the next guess
[42]    →L1 ⍝ Try again
[43]   L2:R←R-1 ⍝ Convert to interest rate
    ∇

Some final notes. These IROR functions have very specific behaviour. You may want to modify them in a number of ways:

1. The iterative process is halted when the new guess (value of i) varies from the previous guess by less than 1E¯11. You may require more or less precision.
2. The process begins with one to three interest rate guesses in the range of 5 to 8 percent. If you have better knowledge of what the unknown interest rate is, you may adjust these initial rates to move them closer to the unknown rate. In this way, you may be able to shorten the process by an iteration or two. In our trials, we found that the process required the same number of iterations regardless of the initial rates … unless the initial rates were very close.
3. In some IROR problems, the interval between cash flows is not uniform. The second cash flow may code in 45 days, the third in 425 days, the fourth in 450 days, and so on. If this is the case, the time period component of the formulas must be changed from the vector ⍳⍴CV to a vector representing fractions of years from the first cash flow to each subsequent cash flow. This change, while straightforward, eliminates the possibility of using cumulative product (×\) in place of the slower exponentiation (*).

2. Limbering Up: Effective Date Searching¶

(The purpose of this column is to work some flab off you APL midsection. Like muscles, your APL skills can atrophy if not exercised with adequate frequency and variety. This column presents a task for you to perform. Set aside a few minutes from your busy schedule and work the task. Mail in your solution and stay tuned for the results.)

Suppose you have two numeric vectors of the same length, IDS and PRICES. IDS is a vector of distinct product ID codes for all available products; PRICES is a vector of corresponding prices for these products. Given a subset list of product ID codes, say LIST, you can easily determine their corresponding prices:

PRICES[IDS⍳LIST]

Suppose prices change over time. Each prices, for a given product, has an associated effective date. For a given product, there may be many different prices, each with a different effective date. To determine the price of a given product, you need to specify not only the product’s ID code, but also the date as of when you want its price.

Here’s an illustration

Your task is to design and define a utility function that takes a list of product Ids and corresponding dates, and returns the prices for those products on those dates. We’re looking for two things: clean design, and program efficiency.

Send your solution to: Vector Production, Brook House, Gilling East, York YO62 4JJ, UK.

The notable functions and their authors’ names will be printed in the next issue of Vector. Good luck and happy limbering.

Ed: When putting this item together I noticed that IROR1[10] would be more efficient if the calculation of present values used decode (⊥) rather than either the exponentiation or the cumulative product techniques. Much to my surprise and horror, the decode algorithm produced a different result to the cumulative product code, but the same result as the exponentiation code. The following demonstrates the problem:

      IROR1 {(⌊¯.05×+/⍵),⍵}?359⍴1000
IROR1[10]
      V+.××⍀(N,2)⍴÷G     ⍝ Cumulative product method
¯53.61702776 ¯4305.479388
      (2 1⍴÷G)⊥⌽V        ⍝ Decode method
¯56.29787915 ¯4736.027327
      (G∘.*-⍳⍴V)+.×V     ⍝ Exponentiation method
¯56.29787915 ¯4736.027327

I can only assume that the cumulative product accumulates errors as the numbers are multiplied together, whereas the other methods do not. In spite of this difference the final IROR result seems to be the same in each case.

Reprinted with kind permission from Zark APL Tutor News, a quarterly publication of Zark Incorporated, 23 Ketchbrook Lane, Ellington, CT06029, USA.