Category: Uncategorized

Update: New Valuation Measure
I wrote about this manuscript in a previous post. It was returned for review after two years. This is a very long time.

I updated it with new data from 2023 and 2024, corrected some misprints, and returned it to editors and reviewers for further review. I hope it gets published soon.

Luckily, this new data did not change any of my conclusions:
- This new valuation measure takes into account the difference between total returns and fundamentals growth.
- This measure works best if we take trailing averaged earnings: Innovations for the autoregressions are independent identically distributed and Gaussian.
- We picked 5 years as averaging window, but this works well also for 10 years, which corresponds to the classic Shiller cyclically adjusted price-earnings ratio (CAPE).
- The classic Shiller CAPE measure shows that the market is currently overvalued: The current ratio is much higher than the historical average or median.
- But our new valuation measure is about equal to the historical average. Thus it shows that the market is fairly valued.
I plan to write yet another manuscript to include annual volatility, see a previous post.
March 14, 2025
S&P 500 Annual Returns with CAPE and Volatility

This is yet another version of the financial simulator. See my GitHub/asarantsev repository CAPE-Volatility.

We have annual volatility and earnings as two factors, as in the previous post. We here use trailing averaged annual earnings over a few years instead of earnings over the last year. This averaging time window ranges from 2 to 10 years. This is similar to the 10-year trailing averaged earnings used in the Shiller cyclically adjusted price-earnings ratio. The inverse of a price-earnings ratio is called earnings yield.

Model: The equations from the previous post stay the same. We model earnings growth normalized by volatility for 1-year earnings. But for the returns, regression has trailing averaged earnings. The data are the same: earnings and inflation for 1927-2024, returns and volatility for 1928-2024. But we add a few data points for earnings and inflation for years 1918-1926. Indeed, we need to compute 10-year averaged earnings, and thus we need a few more trailing data. The data is from my web site.

We take real and nominal earnings and returns. We have 4 types of returns: nominal total, nominal price, real total, real price.

Modification: But then we have a problem with initial conditions. Assume we take a 10-year window. To get earnings yield for year $t$ we need earnings from years $t-9, \ldots, t.$ To get earnings yield for year $t+1$ we need earnings from years $t-8, \ldots, t+1.$ And so on for subsequent years. Thus we need all earnings for each of the last 10 years to be able to simulate returns from this year.

It is not enough to provide only earnings for the last year, or averaged earnings over the last 10 years. These initial conditions would not be sufficient to perform the simulation. We think it would be hard for a user to provide all earnings for each of these 10 years.

Thus I made the following decision. I simply removed the option of choosing initial earnings. The simulator always uses the earnings for the actual last 10 years. The same is true for volatility. We start the simulation from the actual year 2025.

Results for 10-year window: Dependence of returns upon earnings yields is much stronger than in the previous post: $p < 7.5\%$ for all four types of returns. Thus we might state that the 10-year averaging improves the model. However, the $R^2$ values are of the same order: around 20-30%.

Residuals for regressions of returns are normal. Autocorrelation function plots seem better, although still with a strange spike at lag 4.

Simulations show the earnings yield can fluctuate but not nearly as wildly as in the previous post. The picked simulations have yields less than 7%. See the plot below.

Finally, the wealth simulation fluctuates much less. We try 30 years time horizon without any withdrawals or contributions. Overall average returns are more reasonable 9.63%. More importantly, the 90% quantile of average total returns is only 13.65%, much less extreme than in the previous post.

Conclusions: We think the Shiller CAPE, or, equivalently, cyclically adjusted earnings yield with 10 years of averaging window works much better than the annual price-earnings ratio, or equivalently, annual earnings yield. We added annual volatility to the classic Campbell and Shiller’s research and much improved it. Now we have a dynamic stochastic general equilibrium (DSGE) model: We can simulate it for each time step; It is stable in the long run; and it has Gaussian IID innovations.

March 12, 2025

S&P 500 Returns using Earnings Yield & Volatility

In the last few days, I did a lot of research.

In the previous blog post, I have written about the simulator using only one factor: annual volatility, and two equations: for volatility and for total returns. Annual volatility was computed by my undergraduate student Angel Piotrowski as realized standard deviation of daily price returns. I multiplied it by 1000 for normalizing.

Our Model: Here, I extended it to two factors: annual volatility and annual earnings. Previously, annual earnings growth was modeled by another undergraduate student Ian Anderson. We have four equations:

Autoregression of order 1 for logarithmic volatility
Earnings growth divided by volatility as IID
Price returns (excluding dividends) as linear regression with innovations = IID times volatility
Total returns (including dividends) as linear regression with innovations = IID times volatility

These two regressions have factors: this year’s volatility and previous year’s earnings yield: last year’s earnings divided by end of last year’s index. This earnings yield is the classic valuation measure used by financial practitioners.

This model has two versions for the last three equations: nominal (not inflation-adjusted) and real (inflation-adjusted). We use four series of innovations which are multivariate Gaussian: $(G(t), U(t), W(t), Z(t))$ with mean zero. Let us write these equations explicitly. Let $E(t)$ and $V(t)$ be earnings and volatility in year $t.$ Let $S(t)$ be the end of year $t$ index level. Let $Q(t)$ be total returns for this index, including dividends. Then we have:

Volatility: $\ln V(t) = \alpha + \beta \ln V(t-1) + W(t)$

Earnings growth: $\frac1{V(t)}\ln\frac{E(t)}{E(t-1)} = g + G(t)$

Price returns: $\ln\frac{S(t)}{S(t-1)} = a + bV(t) + c\frac{E(t-1)}{S(t-1)} + V(t)U(t)$

Total returns: $Q(t) = A + BV(t) + C\frac{E(t-1)}{S(t-1)} + V(t)Z(t)$

Regression Results: Below is the table. We see that coefficients $a, b, A, B$ are significantly different from zero. In particular, returns and volatility have significant negative dependence. But $c, C$ is not significantly different from zero. Interestingly, we see that dependence of returns upon earnings yield is positive (undervalued markets grow faster) but very weak. A better measure might be cyclically adjusted price-earnings ratio, called Shiller CAPE, see future posts. This ratio is based on 10-year trailing averaged inflation-adjusted earnings.

Presumably this low prediction value is because of high volatility of earnings. For example, in 2008 earnings plummeted so much that the yield plummeted as well, thus making the markets seem overvalued. But in reality, they were undervalued, since earnings rebounded fast, and it took the markets an entire decade to rebound.

Normality: Skewness, kurtosis, and normality tests (Shapiro-Wilk and Jarque-Bera) show that residuals for price and total returns are Gaussian. However, evidence for normality of earnings growth normalized by volatility; or, equivalently, residuals $G(t).$ The Shapiro-Wilk tests for nominal and real earnings growth give us $p$ which are 4.6% and 5.8%. The Jarque-Bera tests give us 0.2% and 0.5%. Kurtosis is especially large: 4.7 for nominal and 4.5 for real versus 3 for Gaussian.

Independence: We show below the autocorrelation functions for original values of residuals: $Z(t)$ and $|Z(t)|.$ We pick the regression for total real returns. Evidence for IID is there but there are problems with this. We are not sure why there is large autocorrelation with lag 4. We did not apply the Ljung-Box omnibus test for the first 5 or 10 autocorrelation values. But we think it would reject the white noise hypothesis. Other regressions for other returns: price real, total nominal, price nominal, show similar autocorrelation plots for residuals and their absolute values.

Next, the autocorrelation function plots for nominal and real growth terms show that these can indeed be modeled by IID. We show plots for real earnings growth. Nominal earnings growth have similar autocorrelation plots.

Simulation. We built a version of the financial simulator which has inputs: initial volatility, initial earnings yield, and time horizon in years. We also allow, as usual, annual contributions or withdrawals. We allow for nominal (not inflation-adjusted) or real (inflation-adjusted) versions. The code and data are available on GitHub/asarantsev repository earnings-yield-annual-simulator. We see the graph of total wealth below. This graph also shows:

No contributions or withdrawals
Real returns
20 years time horizon
Initial volatility: 10.5 (close to historical average)
Initial earnings yield: 7% (close to historical average)

We see that average returns over time and all simulations is 10.21%. This is much higher than historical average returns: 6.7%. This presents a problem: If we start with average data, then we should reproduce historical average returns over many simulations. We also see that with high probability returns are very large: the 90% quantile is 17.11%.

Finally, plot earnings yields for the five chosen simulations. We see that yields can become very high. In two simulations out of five, we have yields greater than 30%. Yields like this do not happen in real life.

What is the reason? We think it is because we have exponential terms in the time series equation for earnings yield $Y(t) = E(t)/S(t).$ Note that the log change in earnings yield can be represented as the difference between price returns and earnings growth:

$\ln Y(t) - \ln Y(t-1) = \ln\frac{E(t)}{E(t-1)} - \ln\frac{S(t)}{S(t-1)}.$

Take the equations for earnings growth and price returns at the beginning of this post. Plug them into the above equation:

$\ln Y(t) - \ln Y(t-1) = -a + (g - b)V(t) - cY(t-1) + V(t)(G(t) - U(t)).$

The change in $\ln Y(t)$ depends on $\ln Y(t-1)$ exponentially. This might lead to volatile fluctuations, especially when the yield is large. Then fluctuations are also large. This needs to be researched further.

Conclusion: We think this makes our research financially unrealistic, despite rigorous statistical analysis. A better way might be to use the CAPE or its inverse, cyclically adjusted earnings yield (using the last few years average for earnings instead of only last year). Another way might be to reproduce my research on the new valuation measure.

March 12, 2025

Financial Simulator, Annual Version

In this post, we discuss a version of this financial simulator for annual data. In the last post, we discussed its currently published monthly version. This is a simplified version based on work by Angel Piotrowski. Here we do not have any bond rates or spreads, any earnings or dividends. We have only annual volatility, previously computed by Angel, and total returns data, available on my web page. The code and data are also available on GitHub/asarantsev repository annual-simulator.

We have log Heston model for annual volatility $V(t)$ :

$\ln V(t) = \alpha + \beta \ln V(t-1) + W(t),\quad \beta = 0.62,\quad \alpha = 0.848.$

Here, we can model $W(t)$ as IID, judging by the ACF plot for $W(t)$ and the ACF plot for $|W(t)|.$

It turns out we cannot quite model $W(t)$ by a normal distribution. But it is quite close. See also the table below. This is quite weak evidence of normality! Next, we model total (including dividends) real (inflation-adjusted) returns $Q(t)$ as follows: $Q(t) = V(t)Z(t)$ for $Z(t)$ IID Gaussian:

These normalized returns $Z(t)$ are very well modeled by IID Gaussian. This is confirmed by the plots above. Similar plots are found for the nominal versions of returns.

We model $(Z(t), W(t))$ as bivariate normal IID: $\mathcal N_2(\mu, \Sigma).$ We found the mean vector and the covariance matrix, for both real and nominal versions. This is available when we run the code.

Next, the code has a function which simulates 1000 paths of wealth:

$\mathcal W(t) = \mathcal W(0)\exp(Q(1)+\ldots +Q(t))$

for $t = 0, \ldots, T.$ This formula would be true for the case when there are no contributions or withdrawals. If each year we have flow (inflow or outflow) $F(t)$ (positive or negative), then the wealth process is

$\mathcal W(t) = \mathcal W(t-1)e^{Q(t)} - F(t)$

The arguments of this function are:

Choose: Real or Nominal?
initial wealth $\mathcal W(0)$
flow $F(t)$ each year
time horizon $T$ in years
initial volatility $V(0)$

We rank 1000 simulations by final wealth. This is the same ranking as by average total returns, if only we have no contributions or withdrawals. But if we do have inflows or outflows every year, then these rankings by wealth vs by returns might be different. For each paths, we compute average total returns $(Q(1) + \ldots + Q(T))/T.$

We pick simulations corresponding to 10%, 30%, 50%, 70%, 90% ranking of final wealth. We chose these because it gives us more or less comprehensive picture of randomness. It is trivial but often neglected that markets have a lot of volatility. It is not enough to provide mean or median returns or wealth. We need the entire distribution.

One could suggest a histogram of final wealth. But a histogram is hard to read by ordinary users. Also, it shows less information than wealth paths. To laypeople, wealth paths fluctuate wildly, and the path which are on top now can drop fast later. But the histogram does not show that.

Why choose these 5 quantiles? We need to show the median (typical wealth) and the outliers. But one should not go too far. Showing the path corresponding to the maximal or minimal final wealth would give a distorted picture: Users might misunderstand these and think that such outcomes are typical. Even showing the paths corresponding to 5% or 95% final wealth are not very typical: These are outliers, atypical outcomes.

The range of 10%, 30%, 50%, 70%, 90% shows full range of outcomes, but only typical outcomes, not outliers.

We also stress that average final wealth does not mean we exceed this wealth in 50% of simulations. The distribution of the final wealth is much skewed, so mean and median are not the same!

Note than in case of withdrawals, the path might end in ruin (bankruptcy; zero wealth). If some of the chosen five paths has this, we show this and we do not compute average total returns. We also compute average total returns for paths which do not end up in ruin.

We compute how many paths end in ruin, thus the empirical probability of ruin. This is a major problem in retirement planning: How much can we withdraw per year so that we do not end up in ruin? The classic withdrawal rate is 4% (of initial wealth) per year, if we adjust these withdrawals for inflation.

We adapted the Python code for running locally, so it contains a function which makes this plot below. This is quite different from the Python code in the simulator, which includes entering data from a web page and putting the output using this web page. We also have HTML code for this web page. For running this Python code locally, obviously we do not need this HTML page, only the Python code itself.

Interestingly, even gigantic 7.5% annual withdrawal rate within 20 years results in only 10.5% ruin probability! Our analysis shows we can be much more permissive with withdrawals than the classic 4% rule.

March 4, 2025

Financial Simulator, Current Version

This is statistical description for the Financial Simulator built using Python Anywhere.

Idea: We divide large and small stock returns by volatility, which makes them closer to normal and independent identically distributed (IID). This is Chicago Board of Exchange Volatility Index (VIX), available from 1986. After division by VIX, large stock returns become IID normal. See the code and data at a GitHub repository.

Data: Try January 1986 – May 2024, $T = 461$ time points, except for volatility, which has $T + 1$ data points. It is taken from Federal Reserve Economic Data (FRED) web site from time series VXOCLS (monthly average data:
Jan 1986 – Feb 1990) and VIXCLS (monthly average data: Mar 1990 – Jun 2024): $V(t),\, t = 1, \ldots, T.$

Large and small stock returns are from Kenneth French Data Library (Portfolios Formed on Size): $Q_0(t),\, t = 1, \ldots, T$ for large stocks (measured by top 30%); and $Q_1(t),\, t = 1, \ldots, T$ for small stocks (measured by middle 40%). These are nominal (not inflation adjusted). To adjust them for inflation, we need to subtract monthly inflation rate, computed from the (not seasonally adjusted) Consumer Price Index taken from FRED. We get then real (inflation-adjusted) versions of returns. We also use short (3-month) and long (10-year) Treasury rates from FRED: end-of-month data Dec 1985 – Jun 2024. We are interested in long-short spread: $S(t), t = 0, \ldots, T.$

Modeling volatility: We model VIX as autoregression of order 1 in the logarithmic scale:

$\ln V(t) = (1 - 0.1181)\ln V(t-1) + 0.3456 + W_0(t),$

where $W_0$ are independent identically distributed innovations with mean zero.
This is a mean-reverting process, stable in the long run. e model residuals $W_0(t)$ as follows: Remove 12 outliers from the 461 in the right tail and model the rest as normal with mean and standard deviation $-0.01593$ and $0.121.$ Thus innovations are modeled as a mixture with weights $p = 12/461, 1-p$ of the uniform distribution upon these 12 outliers, and the above normal distribution. We use kernel density estimation to simulate innovations, with Gaussian kernel and bandwidth derived from the Silverman rule of thumb. In this case the bandwidth is 0.03619. This follows the old blog post.

Modeling long-short bond spread: We model spread as autoregression of order 1 but with additional volatility factor, with residuals normalized by volatility:

$S(t) = - 0.069653 - 0.017407S(t-1) + 0.004575V(t) + V(t)W_1(t)$

This is similar to Ian Anderson’s work but for monthly instead of annual data, and from 1986 instead of 1928.

Modeling stock returns: For nominal returns, we model large and small returns as

$Q_0(t) = 3.583548 - 0.176190 S(t-1) - 0.115088 V(t) + V(t)Z_0(t)$

$Q_1(t) = 3.753727 - 0.038589 S(t-1) - 0.131591 V(t) + V(t)Z_1(t)$

For real returns, we model large and small returns as

$Q_0(t) = 3.581447 - 0.176068 S(t-1) - 0.115111 V(t) + V(t)Z_0(t)$

$Q_1(t) = 3.751626 - 0.038467 S(t-1) - 0.131614 V(t) + V(t)Z_1(t)$

Here $Z_0, Z_1, W_1$ are independent identically distributed trivariate Gaussian innovations,
independent from $W_0$ with mean zero and known covariance matrix. This follows the old blog post, and is similar to Angel Piotrowski’s research. See also large vs small stocks.

March 3, 2025

Earnings Growth/Volatility vs Rates and Spreads

My undergraduate student Ian Anderson did this research available on his GitHub repository. He continued his research on earnings growth from the previous post. In that research, Ian considers growth terms $G(t) = \ln E(t) - \ln E(t-1)$ where $E(t)$ are earnings during year $t.$ These earnings might be nominal (not inflation-adjusted) or real (inflation-adjusted). These are NOT IID Gaussian.

Enter the annual realized volatility computed by another undergraduate student Angel Piotrowski. Denote it by $V(t)$ and divide by it the growth terms. These normalized earnings growth terms $N(t) = G(t)/V(t)$ are, in fact, IID Gaussian.

But are they dependent upon bond rates or spreads? Ian ran simple linear regression of $N(t)$ upon $R(t)$ where we take rate data for January of year $t.$ He has data 1928-2023. All regression residuals have ACF and QQ plots which show they can be modeled by IID Gaussian. Not surprising since terms $N(t)$ are also IID Gaussian. Results for real (inflation-adjusted) earnings are given below.

Quantity	AAA Corporate Rate	10 Year Treasury Rate	1 Year Treasury Rate	AAA – 10YTR Spread
Slope estimate	-1.06	-1.33	-1.36	7.40
Slope p-value	20%	11%	4%	8%
Intercept estimate	10.95	11.33	10.38	-1.46
Intercept p-value	4%	1%	0%	73%
Correlation r	-13%	-17%	-22%	18%

We see that the correlation is significant for 1 year Treasury rate and (to a less extent) for the spread.

It is easy to explain the first correlation: Lower short-term rates make borrowing cheaper and increase access to capital, so earnings grow faster.

What about the second correlation? Such spread is larger in periods of turbulence, with risk premium increasing. But it is not clear why the correlation must be positive.

Results for nominal earnings growth are not very different, except in this case no correlation is strong enough to be statistically significant. The strongest correlation is again with 1-year Treasury rate, which is -14%, and the p-value (the smallest among the four) is 16%.

February 27, 2025

Investment-Grade Corporate Bond Returns 1972-2024

In the two previous blog posts, we considered rates and returns for Bank of America-rated corporate bond portfolios, with ratings AAA, AA, A, BBB, BB, B, CCC. We analyzed annual data and used annual averaged VIX (S&P 500 volatility index) to fit Markov time series models. We were successful: For 5 of these 7 ratings, we created a trivariate model with independent identically distributed trivariate Gaussian innovations. Unfortunately, data was available only starting from 1996. This is less than 30 years.

It would be nice to extend this for a larger time period. We have Moody’s AAA rates back from 1962 (and if you wish only monthly data, all the way from 1919!) But the wealth process (from which we compute total returns) for Bank of America’s corporate bond portfolio is available daily from 1972. To replicate the previous research, we need annual volatility data. But VIX goes only back to 1990, and if we consider volatility for a related index S&P 100, only back to 1986. Thus we need another version of annual volatility. Luckily, Angel Piotrowski provided annual realized volatility for S&P 500 (and its predecessor, S&P 90). This code and data is available on GitHub/asarantsev repository Corporate-Bonds-Annual-Data.

Consider the sequence of rates $R(t)$ which in fact can be modeled by a simple autoregression of order 1:

$R(t) - R(t-1) = a + bR(t-1) + \delta(t).$

Results: $b = -0.06, a = 0.0038$ and the Student T-test gives us $p = 24\%$ . It seems we can model this rate using random walk. The problem of the random walk model is that it is not stable (ergodic) in the long run. The innovations $\delta(t)$ are well modeled as independent identically distributed Gaussian white noise: See the plots below.

This is different from the Bank of America AAA rated bond rates, where we needed to include VIX to make innovations normal.

We also ran autoregression including volatility terms as in previous post. We saw that here, the ACF of absolute values of residuals is not white noise.

But let us continue with total returns computed from wealth process $U(t)$ as follows: $Q(t) = \ln(U(t)/U(t-1)).$ The mean and standard deviation of these terms are $-0.0029$ and $0.071.$

First, we simply test them for independent identically distributed Gaussian: Actually, they are. See the three plots below.

Similarly to Bank of America rated bond portfolio, we write a simple linear regression:

$Q(t) - R(t-1) = k -m(R(t) - R(t-1)) + \varepsilon(t).$

Here $m$ is approximate bond duration. We get $m = 5.8, k = -0.0048.$ This is quite similar to the previous research. But unfortunately, residuals $\varepsilon(t)$ are NOT IID Gaussian. See the three plots below.

We know what to do from our previous research: Rewrite this regression as

$Q(t) - R(t-1) = k - m(R(t) - R(t-1)) + hV(t) + Z(t)V(t).$

Then we get $m = 6.1$ and $k = -0.0076$ and $h = 0.3866.$ Student T-test for $m = 0$ gives very low $p < 0.1\%$ but for $k = 0$ or $h = 0$ we have high p-values, so we fail to reject these null hypotheses. Luckily, these new residuals are, in fact, Gaussian white noise:

The Jarque-Bera and Shapiro-Wilk tests for normality of residuals $\delta(t)$ and $Z(t)$ which give us high p-values (higher than 10%).

I did this work myself, and I am not sure why results here are so different from the previous two blog posts mentioned at the top of this post. This is left for further research.

February 27, 2025

Bank of America Corporate Bonds: Annual Returns

A continuation of research in github.com/asarantsev repository Annual-Bank-of-America-Rated-Bond-Data from my previous post. Consider total returns $Q(t)$ computed from total wealth process $U(t)$ as log change: $Q(t) = \ln U(t) - \ln U(t-1).$ The plot of the wealth process $U(t)$ normalized so that $U(0) = 1$ is given below. We see that high-yield, low-rated bonds provide more long-run returns but with much more risk.

If these were Treasury bonds and there were no risk of default, and if these were zero-coupon bonds (with only principal payment at maturity) then total returns would be equal to the rate minus maturity times rate change. See my manuscript arXiv:2411.03699. The equation is: $Q(t) = R(t-1) - m(R(t) - R(t-1)).$ If the bonds were having coupons, then instead of maturity there would be duration (average time of coupon and principal payments, weighted by payment size). But we add noise (innovation) terms, and an intercept:

$Q(t) - R(t-1) = k - m(R(t) - R(t-1)) + Z(t).$

The maturity is given in the following table, together with analysis of residuals: skewness, kurtosis, Shapiro-Wilk and Jarque-Bera normality test $p$ valued. Also, we take the sum of absolute values of the autocorrelation function for the first five lags, separately for original values of residuals and for their absolute values.

The autocorrelation function plots for $Z(t)$ and for $|Z(t)|$ show that these are independent identically distributed. However, the quantile-quantile plot of $Z(t)$ versus the Gaussian distribution show these are not normal, for most ratings. See the plots below.

This is confirmed by the results of Shapiro-Wilk and Jarque-Bera tests, shown in the table above.

Apply the same technique as in the previous post: Normalize residuals by dividing them by annual average VIX. We get: $Q(t) - R(t-1) = k - m(R(t) - R(t-1)) + V(t)\delta(t).$ We divide this equation by $V(t)$ and get an ordinary least squares regression without intercept. This is not usual, so let us add an intercept:

$Q(t) - R(t-1) = k - m(R(t) - R(t-1)) + hV(t) + V(t)\delta(t).$

Coefficient estimates and analysis of innovations $\delta(t)$ is shown in the table below.

We see the residuals can be well described as Gaussian white noise for ratings BB and higher, especially well for investment-grade bonds. But for B and CCC ratings, not so much. However, judging by the ACF, new residuals (see the second table) are comparable to old residuals (see the first table) in being close to independent identically distributed. See also the following plots for AAA rated bonds:

And for the lowest-rated CCC bonds the situation is different: We see that the first lag is quite significant for both version of the autocorrelation function.

Combining the model above with the results of the previous post, we get the trivariate model:

$\ln V(t) = \alpha + \beta \ln V(t-1) + W(t)$

$R(t) = a + bR(t-1) + cV(t) + V(t)\,\varepsilon(t)$

$Q(t) - R(t-1) = k - m(R(t) - R(t-1)) + hV(t) + V(t)\,\delta(t)$

And the wealth process is given by $U(t) = \exp(Q(1)+\ldots + Q(t)).$ It is possible to show this model for $(V(t), R(t), Q(t))$ is long-term stable and ergodic, because for each of seven ratings, $\beta, b \in (0, 1).$ We have done this for $(V, R)$ in our previous work. For $Q$ this is trivial.

Next, for ratings BB and above, the trivariate innovations sequence $(W, \delta, \varepsilon)$ is modeled as independent identically distributed trivariate Gaussian. Our code find the covariance matrix for these. We do not put it here but an interested reader can run the code.

February 25, 2025

Bank of America Corporate Bonds: Annual Rates

Here we use annual volatility to model bond rates: annual, end-of-year 1996-2024. We take Bank of America bond portfolios with the following seven rates: AAA, AA, A, BBB (investment-grade) and BB, B, CCC (junk, high-yield). Data and code are available on GitHub/asarantsev depository Annual-Bank-of-America-Rated-Bond-Data.

First, consider the rate on the last day of years 1996-2024. Below is the graph of them.

Let $R(t)$ be this rate at end of year $t.$ Model as an autoregression:

$R(t) - R(t-1) = a + bR(t-1) + \delta(t).$

And the results are available in the table below. The last three columns are: Autocorrelation function for innovations, sum of absolute values of the first 5 lags (ACFO); same but for absolute values of innovations (ACFA); Pearson test for $b = 0.$

The autocorrelation function (ACF) plots for $\delta(t)$ and for $|\delta(t)|$ shows that this is well explained by independent identically distributed random variables (white noise). But these are not necessarily normal, judging by the Shapiro-Wilk and Jarque-Bera normality tests. Especially for junk-rated bonds (BB, B, CCC) but also sometimes for investment-grade bonds. The random walk hypothesis could be rejected (using low $p$ values) for all rates (even AAA is just barely above $5\%$ ). See also the plots below. We present only the plots for AAA, other ratings are similar. One can generate these graphs by running the code from the GitHub repository mentioned above.

As usual, we can improve fit and make innovations Gaussian by dividing them by annual volatility. Now we take average annual VIX instead of monthly. This parallels research by Angel Piotrowski mentioned in previous posts. But she computed annual realized volatility, and I use averaged VIX (implied volatility). Let us first fit the log Heston model for VIX 1996-2024:

$\ln V(t) = \alpha + \beta \ln V(t-1) + W(t).$

Here, $\alpha = 1.41$ and $\beta = 1 - 0.475.$ Next, $R^2 = 23.5\%$ and $p = 0.9\%$ for Student $t$ test for $\beta = 1.$ The standard deviation for $W(t)$ is $0.247.$ The normality tests for innovations $W(t)$ give us $p = 27.5\%$ for Shapiro-Wilk and $p = 60.8\%$ for Jarque-Bera. The plots for ACF of $W(t)$ and for $|W(t)|$ show independent identically distributed. See below. Thus the log volatility is indeed modeled by the autoregression of order 1, statistically significantly mean-reverting, with Gaussian innovations. This is similar to Angel Piotrowski’s research.

Consider the autoregression with normalization of innovations $\delta(t)$ by dividing them by volatility $V(t).$ We have then

$R(t) = a + bR(t-1) + V(t)\,\varepsilon(t).$

Divide by VIX and then get an ordinary least squares regression with residuals $\varepsilon(t)$ without intercepts. Let us add intercepts:

$R(t) = a + bR(t-1) + cV(t) + V(t)\,\varepsilon(t).$

Results are available below in the table: Coefficients and analysis of innovations $\varepsilon(t).$

All correlation between $W(t)$ and $\varepsilon(t)$ are not statistically significant. For investment-grade ratings, we have $p$ higher than 5% for all coefficients. But for junk ratings, $p < 5\%$ for $b = 0.$ And for the two bottom ratings, $p < 5\%$ for $c = 0.$ Next, $R^2$ for linear regression with $cV(t)$ for most ratings is much higher than without it. So we need to include this term.

Thus we see a joint model: For $\beta \in (0, 1)$ and $b \in (-1, 0)$ we get:

$\ln V(t) = \alpha + \beta \ln V(t-1) + W(t)$

$R(t) - R(t-1) = a + bR(t-1) + cV(t) + V(t)\varepsilon(t)$

$(W(t), \varepsilon(t)) \sim \mathcal N_2([0, 0], \Sigma)$ IID

As discussed, we might consider $\Sigma$ to be the diagonal matrix, but might as well make it a complete matrix. This model fits very well. A disadvantage is that we have only ~30 years of data. We do need to normalize innovations of rates by dividing these by VIX. We do need the term $cV(t).$

In our previous research, we proved long-term stability of this bivariate model. This is true not just for Gaussian innovations, but for more general cases, under certain conditions.

We continue this research in the next post, where we model total returns.

February 25, 2025

New Bubble Measure Work Replicated

This is the work with my undergraduate student Angel Piotrowski. She used her annual volatility data 1928-2023 to replicate my previous work published on arXiv and discussed in my previous blog entry about the new valuation measure for the Standard & Poor 500 and its predecessors. She used the end-of-year S&P 500 close trading price instead of the January S&P 500 daily close averaged price used by me. She used December Consumer Price Index (CPI) data instead of me using January CPI. And she used data only 1928-2023 instead of my work 1871-2023. But this almost a century of data is still enough to make conclusions.

First, adjust earnings for inflation and consider the trailing average $E(t)$ of earnings for the last 5 years. This is similar to the classic Campbell-Shiller approach, when the cyclically adjusted price-earnings ratio (CAPE) features last 10 years of earnings. But we chose 5 years to make room for more data. Next, take inflation-adjusted wealth $S(t)$ at end of year $t.$ Consider the linear regression

$\ln S(t) - \ln E(t) - ct = a + b(\ln S(t-1) - \ln E(t-1) - c(t-1)) + \delta(t).$

This has the meaning of subtracting the trend from relative growth of wealth over earnings. Historical earnings averaged 1-2% per year and wealth growth is around 6-7% per year so the value for the trend $c$ must be 4-5%. This equation shows that after detrending, this is the classic autoregression of order 1. We can rewrite this in the more standard ordinary least squares regression form:

$Q(t) = a + bc + t(c - bc) + (b-1)(\ln S(t-1) - \ln E(t-1)) + \delta(t),$

where we define for short notation the quantity which we called implied dividend yield:

$Q(t) := \ln S(t) - \ln E(t) - \ln S(t-1) - \ln E(t-1)$

This gives us $a + bc = 0.1351$ and $c - bc = 0.0101$ and $b - 1 = -0.1307.$ All three coefficients are significantly different from zero: Student T-test gives $p$ values $0.1\%, 2.9\%, 1.9\%.$ But the Jarque-Bera test shows that innovations are not Gaussian: $p = 0.1\%.$ This is confirmed by the following quantile-quantile plot below. But the autocorrelation function plots for $\delta(t)$ and for $|\delta(t)|$ below show that these are independent identically distributed.

This stays in contrast with our original research, when residuals (innovations) are independent identically distributed and Gaussian. Further research will include extending this to years 2024 and 1924-1928 when we have total returns and volatility data for shifted years.

The GitHub repository New Valuation Measure Replication contains code, data, and these graphs.

February 20, 2025

	Point Estimate	95% Confidence Interval	p-value for Student test
Real Returns $a$	0.1655	[0.064, 0.267]	0.002
Real Returns $b$	-0.0141	[-0.024, -0.005]	0.004
Real Returns $c$	0.1453	[-0.903, 1.194]	0.784
Real Returns $A$	0.1668	[0.068, 0.265]	0.001
Real Returns $B$	-0.0136	[-0.023, -0.004]	0.004
Real Returns $C$	0.5494	[-0.465, 1.564]	0.285
Nominal Returns $a$	0.1693	[0.074, 0.264]	0.001
Nominal Returns $b$	-0.0133	[-0.022, -0.004]	0.003
Nominal Returns $c$	0.4432	[-0.537, 1.424]	0.372
Nominal Returns $A$	0.1706	[0.079, 0.262]	0.000
Nominal Returns $B$	-0.0127	[-0.021, -0.004]	0.004
Nominal Returns $C$	0.8473	[-0.097, 1.792]	0.078

Data	Skewness (normal =0)	Kurtosis (normal = 0)	Shapiro-Wilk $p$	Jarque-Bera $p$
AR(1) Innovations for Log Volatility	0.59	0.057	0.94%	6.1%
Normalized Total Nominal Returns	0.22	-0.21	11%	60%
Normalized Total Real Returns	0.24	0.023	20%	63%

Rating	$m$	Skewness	Kurtosis	Shapiro-Wilk $p$	Jarque-Bera $p$	ACF of $Z(t)$	ACF of $\|Z(t)\|$
AAA	6.03	-2.049	5.173	0.014%	<0.001%	0.539	0.687
AA	4.89	-0.941	1.154	2.549%	5.831%	0.944	0.874
A	4.94	-0.894	1.307	7.877%	5.726%	0.977	0.704
BBB	5.17	-0.573	0.027	25%	46%	0.38	0.604
BB	3.81	-1.75	3.64	0.054%	<0.001%	0.83	0.245
B	3.12	-2.036	4.323	0.003%	<0.001%	1.17	0.653
CCC	2.55	-2.304	5.34	0.001%	<0.001%	0.9	0.716

Rating	$k$	$m$	$h$	Skewness	Kurtosis	Shapiro-Wilk $p$	Jarque-Bera $p$	ACF of $Z(t)$	ACF of $\|Z(t)\|$
AAA	0.0661	7.0787	-0.0034	-0.737	0.576	28%	23%	0.922	0.361
AA	0.0453	5.3226	-0.0023	0.126	-0.039	46%	96%	0.57	0.761
A	0.0423	5.3423	-0.0022	-0.181	0.02	35%	93%	0.73	0.459
BBB	0.0293	5.6074	-0.0016	-0.232	-0.775	67%	62%	0.875	0.498
BB	0.0422	3.6671	-0.0024	-0.894	1.426	19%	6.3%	0.574	0.888
B	0.0682	2.9970	-0.0050	-1.562	3.479	0.397%	<0.001%	0.998	0.354
CCC	0.0712	2.6532	-0.0075	-1.588	2.649	0.072%	0.005%	0.835	0.723

Rate	$b$	$a$	Stdev of residuals	Skew	Kurtosis	Shapiro-Wilk $p$	Jarque-Bera $p$	ACFO	ACFA	Pearson Test $p$
AAA	-0.21	0.0081	0.008	0.994	1.241	0.029	0.041	0.629	0.627	0.051
AA	-0.23	0.0087	0.009	0.747	0.842	0.108	0.18	0.716	0.918	0.049
A	-0.26	0.011	0.01	0.571	0.53	0.138	0.397	0.48	0.744	0.043
BBB	-0.33	0.017	0.012	0.799	1.671	0.09	0.044	0.262	0.699	0.025
BB	-0.46	0.031	0.02	1.479	3.951	0.007	<0.001	0.482	0.814	0.009
B	-0.57	0.048	0.026	1.612	3.629	0.003	<0.001	0.366	0.727	0.003
CCC	-0.57	0.082	0.055	1.338	2.051	0.007	0.001	0.576	0.626	0.004

Rate	$a$	$b$	$10000*c$	Stdev	Skew	Kurt	S-W	J-B	ACFO	ACFA
AAA	0.01	-0.14	-2.41	0.00037	0.774	0.664	22%	19%	0.352	0.816
AA	0.0099	-0.116	-2.78	0.0004	0.67	0.807	43%	24%	0.402	0.924
A	0.0087	-0.157	-1.06	0.00044	0.342	0.2	93%	74%	0.395	0.963
BBB	0.0099	-0.279	-2.1	0.00051	0.039	-0.4	96%	91%	0.405	0.779
BB	0.015	-0.539	10.0	0.00076	0.303	-0.21	88%	79%	0.656	0.594
B	0.0235	-0.77	20.3	0.00096	0.35	-0.126	46%	74%	0.576	0.722
CCC	0.0217	-0.73	41.5	0.00198	0.58	0.25	36%	44%	0.817	0.948