Replicate average risk-free rate in Wachter 2005, "Solving models with external habit" I might be making a really simple mistake somewhere, but I thought I'd ask anyway. I'm trying to replicate the results in Wachter 2005, &quot;Solving models with external habit&quot;. (You can also find the prepublication version here for free.) However, I'm getting a slightly different result for some quantities and I'm wondering if I'm interpreting the parameters correctly. In particular, when I input the parameters listed in Table 1, I get a different average risk-free rate in both the &quot;CC value&quot; case and the &quot;Wachter value&quot; case.

From Tables 2 and 3, we see that the &quot;CC value&quot; case should give an average risk-free rate, $E[r^f]$, of $0.94\%$ and the
&quot;Wachter value&quot; case should give $1.47\%$.
The value of the average risk-free rate can be derived analytically, and is given in equation (11):
$$
E[r^f_t] = - \log \delta + \gamma g - \frac{\gamma (1-\phi) - b}{2}. \tag{11a}
$$
When I substitute in the parameter values from Table 1, I get
$0.47\%$ for the CC parameterization and
$2.47\%$ for the Wachter parameterization.
I wonder if I am annualizing correctly? I'm somewhat confident that I'm inputting the parameters correctly, because the paper does show that for the Campbell Cochrane calibration (CC values), a derived parameter called $S_{\text{max}}$ should be equal to $0.0939$,
where
\begin{align}
\bar S &amp;= \sigma_\nu \sqrt{\frac{\gamma}{1-\phi - b/\gamma}} \tag{5}\\
s_{\text{max}} &amp;= \bar s + \frac 12 (1 - \bar S^2), \tag{6}
\end{align}
Lower case letters denote the natural log of the capitalized variable. When I substitute in the &quot;CC values&quot;, I get $S_{\text{max}} = 0.09384$. So, this seems correct. The paper implicitely reports this parameter in a footnote:

Any thoughts as to why I'm getting different values for the average risk-free rate?
I've included Python code below to show what I'm doing.
# Campbell Cochrane parameters
ann_fact = 12 #annualization factor. Monthly
g = 1.89 / ann_fact / 100
sigma_v = 1.50 / np.sqrt(ann_fact) /100
gamma = 2.00
b = 0.00
phi = (0.87)**(1/ann_fact)
delta = (0.90)**(1/ann_fact)

# ## Wachter 2005 Parameters
# ann_fact = 4 #annualization factor. Quarterly
# parameterization_name = 'Wachter (2005)'
# ann_fact = ann_fact
# g = 2.20 / ann_fact / 100
# sigma_v = 0.86 / np.sqrt(ann_fact) / 100
# gamma = 2.00
# b = 0.011
# phi = (0.89)**(1/ann_fact)
# delta = (0.93)**(1/ann_fact)

bar_S = sigma_v * np.sqrt(gamma/(1 - phi - b/gamma))
bar_s = np.log(bar_S)
s_max = bar_s + 1/2 * (1 - bar_S**2)
Smax = np.exp(s_max) 
print(&quot;Smax &quot;, Smax) # CC1999 should be 0.0939
r_f = (-np.log(delta) + gamma * g - 1/2 * (gamma * (1-phi) - b))
print(&quot;Average risk-free rate&quot;, r_f * ann_fact * 100)
# CC1999 should be 0.94 %
# Wachter 2005 should be 1.47 %
I think the answer just comes down to rounding error. The reported values of $\delta$ are annualized. In the case of the CC1999 calibration, the annualized $\delta$ is $0.90$ and the implied monthly is $\delta = 0.991258$. Because of the compounding over 12 periods (4 in the Wachter 2005/2006 case), the average risk-free rate is very sensitive to small changes in this value. Taking the other parameters as given and solving (11a) for the appropriate $\delta$ in each calibration leads to the following annualized values that are needed to reproduce the desired risk-free rates:
\begin{align}
\delta &amp;= 0.8957829809431154 \text{ (instead of 0.90 in CC calibration)} \\
\delta &amp;= 0.9384246954559854 \text{ (instead of 0.93 in Wachter calibration)}
\end{align}
Below is the code to verify this.
# Campbell Cochrane parameters
ann_fact = 12 #annualization factor. Monthly
g = 1.89 / ann_fact / 100
sigma_v = 1.50 / np.sqrt(ann_fact) /100
gamma = 2.00
b = 0.00
phi = (0.87)**(1/ann_fact)
# delta = (0.90)**(1/ann_fact)
delta_CC = (np.exp(- 0.94/ann_fact/100 + (gamma * g - 1/2 * (gamma * (1-phi) - b)) ))**ann_fact
delta = (delta_CC)**(1/ann_fact)

# ## Wachter 2005 Parameters
# ann_fact = 4 #annualization factor. Quarterly
# parameterization_name = 'Wachter (2005)'
# ann_fact = ann_fact
# g = 2.20 / ann_fact / 100
# sigma_v = 0.86 / np.sqrt(ann_fact) / 100
# gamma = 2.00
# b = 0.011
# phi = (0.89)**(1/ann_fact)
# # delta = (0.93)**(1/ann_fact)
# delta_wachter = (np.exp(- 1.47/ann_fact/100 + (gamma * g - 1/2 * (gamma * (1-phi) - b)) ))**ann_fact
# delta = (delta_wachter)**(1/ann_fact)

bar_S = sigma_v * np.sqrt(gamma/(1 - phi - b/gamma))
bar_s = np.log(bar_S)
s_max = bar_s + 1/2 * (1 - bar_S**2)
Smax = np.exp(s_max)
print(&quot;Smax &quot;, Smax) # CC1999 should be 0.0939
r_f = (-np.log(delta) + gamma * g - 1/2 * (gamma * (1-phi) - b))
print(&quot;Average risk-free rate&quot;, r_f * ann_fact * 100)
# CC1999 should be 0.94 %
# Wachter 2005 should be 1.47 %
```