Issue with the mm_root function in computing the implied cost of equity capital

[Paul Demere]

Hello Stata Users,
I have encountered an issue where the mm_root command (part of the moremata package) fails to return anything when I use upper and lower bounds of 0 and 1, even though (a) these are the bounds that are theoretically required and (b) using Excel solver for a few observations shows that the solution for the implied cost of equity capital falls within these bounds. Further, when expanding the upper and lower bounds out (I have tried between -10 and 10, as well as -30,000 to 30,000 and a variety of bounds in-between) I will get implied cost of equity capital estimates returned, but every single observation is either at the lower or the upper bound. Varying the tolerance from 0 to 1 and increasing the number of iterations to 10,000 does not help. Is there something I am missing, or is there another better way to code this problem?

For reference, I am trying to compute the cost of equity capital in accordance with Gebhardt, Lee, and Swaminathan (2001), as modifed by Li and Mohanram (2014) (see below for complete citations). The specific code I am using is:

mata
z=J(1,1,.)
st_view(z,., "glsicc p3prcc bve fore1 fore2 fore3 fore4 fore5 indmedroe dvrate")
function w(x,p,b0,f1,f2,f3,f4,f5,ii,dv) {
    tf = f1+f2+f3+f4
    fi = f5-ii
    return(-p+b0+(f1-x*b0)/(1+x)+(f2-x*(b0+f1*dv))/(1+x)^2+(f3-x*(b0+dv*(f1+f2)))/(1+x)^3+(f4-x*(b0+dv*(f1+f2+f3)))/(1+x)^4+(f5-x*(b0+dv*tf))/(1+x)^5+((f5-fi/7)-x*(b0+dv*(tf+f5)))/(1+x)^6+((f5-fi*2/7)-x*(b0+dv*(tf+2*f5-fi/7)))/(1+x)^7+((f5-fi*3/7)-x*(b0+dv*(tf+3*f5-fi*3/7)))/(1+x)^8+((f5-fi*4/7)-x*(b0+dv*(tf+4*f5-fi*6/7)))/(1+x)^9+((f5-fi*5/7)-x*(b0+dv*(tf+5*f5-fi*10/7)))/(1+x)^10+((f5-fi*6/7)-x*(b0+dv*(tf+6*f5-fi*15/7)))/(1+x)^11+(ii-x*(b0+dv*(tf+7*f5-fi*3)))/((1+x)^12)+(ii-x*(b0+dv*(tf+7*f5-fi*3)))/(x*(1+x)^12))
}
for (i=1;i<=rows(z);i++) {
    r=mm_root(glsicc=.,&w(),0,1,1e-9,1000,z[i,2],z[i,3],z[i,4],z[i,5],z[i,6],z[i,7],z[i,8],z[i,9],z[i,10])
    z[i,1]=glsicc
}
end

References:
Gebhardt, W., C. Lee, and B. Swaminathan. 2001. "Toward an Implied Cost of Capital". Journal of Accounting Research 39(1): 135-176.
Li, K. and P. Mohanram. 2014. "Evaluating Cross-Sectional Forecasting Models for Implied Cost of Capital" Review of Accounting Studies 19:1152-1185.

[Steve Samuels]

This code is complicated already, but the long lines makes it difficult to parse. Please use

#delim;

and respost with shorter lines.

[Hua Peng (StataCorp)]

It is impossible to give help without the dataset for parameters "glsicc p3prcc bve fore1 fore2 fore3 fore4 fore5 indmedroe dvrate". Even with the dataset, it will be difficult since -moremata- is a compiled Mata library and users do not have access to the source code. I would suggest you contact Ben Jann, the author of -moremata-, directly.

[Christophe Kolodziejczyk]

You can look at the source code of the moremata package by typing help moremata_source. It is actually very instructive on Mata generally.

That being said you have very complicated function so it could be that the problem comes from here. An error in the coding of the function could mean that the problem has no solution. I would try first with a simpler function and with a dataset delivered with Stata. If you have problems there send the code and a log file of what has gone wrong.

[Hua Peng (StataCorp)]

Thanks Christophe Kolodziejczyk, I did not know -help moremata_source-.

[Paul Demere]

Thank you very much for the replies! On the advice of Hua Peng (StataCorp), I reached out to Ben Jann. He pointed out that the issue was with the bounds I was using, as the function is undefined at the bounds. By replacing the lower bound with smallestdouble() and the upper bound with 1-epsilon(1), I was able to run the code successfully. Included below is the code modified with the new bounds and made more readable at the request of Steve Samuels.

generate glsicc=.
mata
z=J(1,1,.)
st_view(z,., "glsicc p3prcc bve fore1 fore2 fore3 fore4 fore5 indmedroe dvrate")
function w(x,p,b0,f1,f2,f3,f4,f5,ii,dv) {
    tf = f1+f2+f3+f4
    fi = f5-ii
    t1 = (f1-x*b0)/(1+x)
    t2 = (f2-x*(b0+f1*dv))/(1+x)^2
    t3 = (f3-x*(b0+dv*(f1+f2)))/(1+x)^3
    t4 = (f4-x*(b0+dv*(f1+f2+f3)))/(1+x)^4
    t5 = (f5-x*(b0+dv*tf))/(1+x)^5
    t6 = ((f5-fi/7)-x*(b0+dv*(tf+f5)))/(1+x)^6
    t7 = ((f5-fi*2/7)-x*(b0+dv*(tf+2*f5-fi/7)))/(1+x)^7
    t8 = ((f5-fi*3/7)-x*(b0+dv*(tf+3*f5-fi*3/7)))/(1+x)^8
    t9 = ((f5-fi*4/7)-x*(b0+dv*(tf+4*f5-fi*6/7)))/(1+x)^9
    t10 = ((f5-fi*5/7)-x*(b0+dv*(tf+5*f5-fi*10/7)))/(1+x)^10
    t11 = ((f5-fi*6/7)-x*(b0+dv*(tf+6*f5-fi*15/7)))/(1+x)^11
    t12 = (ii-x*(b0+dv*(tf+7*f5-fi*3)))/((1+x)^12)
    tv = (ii-x*(b0+dv*(tf+7*f5-fi*3)))/(x*(1+x)^12)
    return(-p+b0+t1+t2+t3+t4+t5+t6+t7+t8+t9+t10+t11+t12+tv)
}
for (i=1;i<=rows(z);i++) {
    r=mm_root(glsicc=.,&w(),smallestdouble(),1-epsilon(1),1e-9,1000,z[i,2],z[i,3],z[i,4],z[i,5],z[i,6],z[i,7],z[i,8],z[i,9],z[i,10])
    z[i,1]=glsicc
}
end

[Julia Schmidt]

Dear STATA Users,
I have a similar problem as Paul. I tried to estimate the ICC defined by Claus and Thomas (2001) as implemented by Li and Mohanram (2014). I used a similar approach as in Paul’s first posting and my results are mostly 0 or 1. There are some exceptions where STATA generates the correct value.
So I tried to change my STATA code according to your solution, Paul, but still, the results are not getting better.
Surprisingly, wenn I estimate ICC_GLS with your code, Paul, it generates the correct proxies.
I hope, someone can help me and finds the mistake in my code:

generate icc_ct=.
mata
z=J(1,1,.)
st_view(z,., "icc_ct price bv_t0 feps_t1 feps_t2 feps_t3 feps_t4 feps_t5 growth payoutratio")
function y(x,p,b0,f1,f2,f3,f4,f5,g,dv) {
    tf = f1+f2+f3+f4
    t1 = (f1-x*b0)/(1+x)
    t2 = (f2-x*(b0+f1*(1-dv)))/(1+x)^2
    t3 = (f3-x*(b0+(1-dv)*(f1+f2)))/(1+x)^3
    t4 = (f4-x*(b0+(1-dv)*(f1+f2+f3)))/(1+x)^4
    t5 = (f5-x*(b0+(1-dv)*tf))/(1+x)^5
    tv = ((f5-x*(b0+(1-dv)*tf))*(1+g)/((x-g)*(1+x)^5)
    return(-p+b0+t1+t2+t3+t4+t5+tv)
}
for (i=1;i<=rows(z);i++) {
    r=mm_root(icc_ct=.,&y(),smallestdouble(),1-epsilon(1),1e-9,1000,z[i,2],z[i,3],z[i,4],z[i,5],z[i,6],z[i,7],z[i,8],z[i,9],z[i,10])
    z[i,1]=icc_ct
}
end

Thanks a lot!
Julia


References:

• Claus, James and Jacob Thomas (2001): Equity Premia as Low as Three Percent? Evidence from Analysts’ Earnings Forecasts for Domestic and International Stock Markets, in: The Journal of Finance, Vol. 56, pp. 1629 – 1666.
• Li, Kevin K. and Partha Mohanram (2014): Evaluating Cross-Sectional Forecasting Models for Implied Cost of Capital, in: Review of Accounting Studies, Vol. 19, pp. 1152 – 1185.

[Ingrid Qiu]

Hi Stata Users,
I encountered a similar question when using mm_root(). I tried the method suggested above but it doesn't work for me. But after I changed the initial value, mm_root() returns the correct answers. I feel like mm_root() is very sensitive to the bound and the initial values. I am wondering whether there are something we need to pay attention to when we use this function. Thanks!

[Hua Peng (StataCorp)]

Julia Schmidt, your function is different from Paul's, it looks like the expansion of the same function but yours stops at 5th or 6th order, Paul's goes to 12th or 13th order, depending what the underlying function behaves given different coefficients "icc_ct price bv_t0 feps_t1 feps_t2 feps_t3 feps_t4 feps_t5 growth payoutratio", they may have very different roots. Again, this is just speculation without looking the coefficients.

Ingrid Qiu, I would be very interested if you could post the code of the problem you encountered.

[Julia Schmidt]

Hua Peng (StataCorp) The function is different because I try to estimate the ICC following Claus and Thomas, while Paul follows Gebhardt et al. The approaches differ in terms of the periods considered and how the terminal value is calculated. Since he was describing similar issues in his first post, I was hoping that someone could tell me, how to change my function to get the real results and not just the boundaries.
Ingrid Qiu what exactly is the initial value and how did you change it?
Thanks for your help

[Hua Peng (StataCorp)]

In your function, the last term is

tv = ((f5-x*(b0+(1-dv)*tf))*(1+g)/((x-g)*(1+x)^5)

which means the function is not defined when x = g, i.e., at the value of the growth variable. In Paul's version, the function is not defined at x=0. Hence you might need try the bounds from smallestdouble() to z[i,9]-epsilon(1); and z[i,9]+epsilon(1) to 1 if the value of growth variable is between 0 and 1.

[Ingrid Qiu]

Hua Peng (StataCorp) Julia Schmidt Thank you for your interest.

use data

mata:

st_view(Y=.,.,"Y")
st_view(val=.,.,"val")

P = 0.9676741647
t = 0.5
function eq1(y,Y,val,P,t) return(sum((Y:<y):*val)/sum(val) - (t-1+P)/P)
rc = mm_root(x1=.,&eq1(),0,9,4.5,1000,Y,val,P,t)

function eq2(y,Y,val,P,t) return(sum((Y:<y):*val)/sum(val) - t/P)
rc = mm_root(x2=.,&eq2(),0,9,0,1000,Y,val,P,t)

The above code gives the correct solution. I used

  quant = mm_quantile(Y,1,t)

as the initial value, but this will yield a wrong solution. Also, if I change the interval to 6-9, mm_root() will also return the wrong solution. I am not sure what makes this problem.

[Hua Peng (StataCorp)]

The problem is the functions eq1() and eq2() are not continues in the domain you specified. Note Y:<y is discontinues at y. The following shows eq1() changes around 0.02 when the input is changed of 1e-15 at Y[38, 1]. Actually, function eq1() is discontinues at any observation Y where variable val is not 0.

. mata:
------------------------------------------------- mata (type end to exit) -----
: x = Y[38, 1]

: x
  6.420057774

: r1 = eq1(x, Y, val, P, t)

: r1
  -.2658067825

:
: r2 = eq1(x+(1e-15), Y, val, P, t)

: r2
  -.2411065604

: end
-------------------------------------------------------------------------------

mm_root() implements Brent's method https://en.wikipedia.org/wiki/Brent's_method which does assume the underlying function to be continues. That said, the algorithm is robust enough and can often deal with discontinues functions which behave reasonably. But the results need be carefully examined. In your above code, the first call to mm_root() is using a tolerance of 4.5 in a range 0-9 which seems way too large. Change it to 1e-5 gives the following:

: rc = mm_root(x1=.,&eq1(),0,9,1e-5,1000,Y,val,P,t)

: x1
  6.645090348

: rc
  0

:
: //Plot it back
: y = eq1(x1, Y, val, P, t)

: y
  -.0004896063

Use the same parameters but change the search interval to 6-9 give the same result:

: rc = mm_root(x1=.,&eq1(),6,9,1e-5,1000,Y,val,P,t)

: x1
  6.645088397

: rc
  0

:
: //Plot it back
: y = eq1(x1, Y, val, P, t)

: y
  -.0004896063

mm_root() believes it finds the root since rc is 0 and the result seems to reasonable when plotting back.

For eq2(),

: rc = mm_root(x1=.,&eq2(),0,9,1e-5,1000,Y,val,P,t)

: x1
  6.657811668

: rc
  0

:
: //Plot it back
: y = eq2(x1, Y, val, P, t)

: y
  -.0102072419

Lowering the tolerance and increasing the iteration do not improve the result:

: rc = mm_root(x1=.,&eq2(),0,9,1e-15,10000,Y,val,P,t)

: x1
  6.657814026

: rc
  0

:
: //Plot it back
: y = eq2(x1, Y, val, P, t)

: y
  -.0102072419

mm_root() believes it finds the root since rc is 0, but the results from plotting back raises concern which may be caused by the discontinuity of the underlying function.

[Julia Schmidt]

Hua Peng (StataCorp) Thank you so much!! Now it's working

[Attaullah Shah]

@ Julia Schmidt Can you please share an example data from your data and the code that worked for you.