Agent-based model generating stylized facts of fixed income markets

Abstract

We develop an agent-based model (ABM) of a financial market with multiple assets belonging either to the fixed income or equity asset classes. The aim is to reproduce the main stylized facts of fixed income markets with regards to the emerging dynamics of the yield curves. Our ABM is rooted in the market model of Kaizoji et al. (J Econ Behav Organ 112:289–310, 2015) formulated with two types of traders: the rational and risk-averse fundamentalist investors and the noise traders who invest under the influence of social imitation and price momentum. The investors involved in the present market model diversify their investments between a preferred stock equivalent to a perpetual bond and multiple bonds of selected maturities. Among those, a zero-coupon bond provides a constant rate of return, while the prices of the coupon-paying bonds are determined at each time step by the equilibrium between the investors’ demands and supplies. As a result, the ABM creates an evolving yield curve determined by the aggregate impact of the traders’ investments. In agreement with real markets, it also produces transient turbulent periods in the prices’ time series as well as a humped term structure of volatility. We compare the dynamics arising from different processes governing the risk-free rate with those of the historical US Treasury market. Introducing Vasicek’s model of interest rates to both synthetic and empirical rates demonstrates the capacity of our ABM in reproducing the main characteristics of the surface of autocorrelation of the volatilities of the yields to maturity of the US Treasury bonds for the selected time-frame.

1 Introduction

This contribution aims at unveiling the potential residing in the creation of agent-based models (ABMs) of the fixed income asset class involving fundamentalist and chartist traders. To the best of the authors’ knowledge, we provide the first ABM of a financial market containing multiple assets belonging to the fixed income asset class. The value of this type of computational economic model lies in its capacity to produce transient market phenomena originated by the interactions of the individual agents. As pointed out in Sornette (2014), these models are freed from any stationary condition and can reproduce such periods of market unrest instigated by prices’ rallies and their subsequent corrections. Examples of ABMs developed to describe stock markets can be found in Kyle (1985), Black (1986) and Samanidou et al. (2007), among others.

The literature provides a significant number of equilibrium models implemented to analyse the dynamics of yield curves and to forecast their evolution through the behaviour of the corresponding forward rates. Such analyses are, e.g. found in Nelson and Siegel (1987), Cox et al. (1985), Vasicek (1977), Duffie and Kan (1996), Dai and Singleton (2002), Duffee (2002) and Heath et al. (1992). These models have, however, proven limited with regards to their capacity to account for out-of equilibrium conditions corresponding to the ever-recurring extreme transient phenomena observed in reality. Moreover, they do not consider the relationship between the equity and fixed income asset classes. The market model developed below aims at laying a foundation towards closing this gap. We thus propose a new approach to generate stylized facts of fixed income markets by a market model including the two previously mentioned asset classes. As exposed in Kaldor (1961), stylized facts represent “broad tendencies” understood as empirical truths observed in the markets of interest. The focus is therefore drawn on the emerging dynamics associated with the yields to maturity and evolving as a function of the traders’ investments updated at each time step in conformity with the rules developed in this novel type of ABM.

The present formulation derives from the market model proposed by Kaizoji et al. (2015) introducing two types of agents investing in two assets. The risk-free asset provides a constant rate of return and the remaining asset’s price is subjected to the market clearing process applied at each time step. This risky asset also pays dividends determined from a stochastic multiplicative process. The agents are either fundamentalists or noise traders and do not change their strategies in time. This is in contrast to earlier ABMs in which the traders may switch between predictor (Arthur et al. 1996) or between chartist and fundamentalist strategies (De Grauwe et al. 1995; Brock and Hommes 1997). The fundamentalist strategy is as described in Brock and Hommes (1998) and in Chiarella et al. (2009). The chartist strategy is derived from the developments of Lux and Marchesi (1999). In a nutshell, their investments are updated probabilistically at each time step under the influence of other agents and of the momenta in the prices’ time series. At the aggregated level, Kaizoji et al. (2015)’s model proved able to reproduce fat-tail distributions of returns, slow decaying autocorrelations of absolute returns, fast decaying autocorrelations of signed returns, volatility clustering and transient faster-than-exponentially growing prices associated with bubbles. A market model characterized by the interactions of fundamentalist and chartist traders through the price-vector of dimension higher than two is provided in Xu et al. (2014). Other models introducing heterogeneous agents investing in a multi-asset market can be found in Borghesi and Bouchaud (2007), Chiarella et al. (2007), Fedyk et al. (2013) and Eckrot et al. (2016). More recent articles revealing a renewed interest in agent-based modelling are proposed by Gualdi et al. (2015), Baghestanian et al. (2015), Bouchaud (2013) and Hommes and LeBaron (2018). Concerning the application to fixed income, Braun-Munzinger et al. (2018) in particular present a model involving three distinct agents investing in a unique corporate bond.

Section 2 introduces the market framework involving the fixed income and equity asset classes as well as the two trader types investing in this framework. The market clearing process is defined subsequently before the unveiling of the parameter selection and the initialization of the variables achieved in Sect. 3. The latter section concludes with the analyses of the time series generated in two representative simulations, the first one corresponding to a quiet market regime and the other one exhibiting the characteristics of a turbulent regime. Section 4 presents the analyses comparing the dynamics emerging from the market model with those of the US Treasury market realized between November 1993 and June 2020. The section starts with the implementation of Vasicek’s model of interest rates to calibrate the synthetic risk-free rates of the simulations. It continues with the presentation of the term structures of volatility and surfaces of autocorrelations of the yields to maturity associated with each model set-up and with the US Treasury market for comparison. Section 5 concludes.

2 Model set-up

Our market model is made of a stock and of multiple bonds of selected maturities. One of the bonds is a zero-coupon bond providing a constant reference rate throughout the simulations. The agents belong either to the fundamentalist or to the noise trader type. The fundamentalist formulates expectations of the future returns as well as of the risk associated with each asset. The fundamentalist’s investments are ruled by the need to maximize a constant relative risk aversion (CRRA) expected utility function. The noise traders’ investments are subjected to social imitation, adherence to momentum trading and to the influence of a time-varying herding propensity parameter. Their individual investments are updated probabilistically at each time step. Each trader type is represented by a unique agent expressing the aggregate excess demands for each asset at each time step. The market clearing condition imposing the equilibrium between the supplies and demands for the assets determines the endogenous asset prices.

In the following, we start with presenting the market framework considering the fixed income asset class. The assets’ and wealth dynamics are then derived. The description of both aggregated traders is then given, starting with the fundamentalist and the associated generalized optimization problem. The threefold investment process governing the chartist’s investments is subsequently defined before the presentation of the market clearing achieved by the Walrasian auctioneer (Walras 1874). The latter process provides the system of nonlinear equations ruling the updates of the asset prices at each time step. The selection of the parameter values initializing the market model follows, before the presentation of the time series generated by a configuration of the market model considering six coupon paying bonds and a unique preferred stock.

2.1 The fixed income market framework

Our aim is to constrain the model to involve traders following only “aggressive” investment strategies. The selection of such traders is justified by their catalysing effect on the market through the permanent rebalancing of their portfolio. This is in contrast to passive investors who “freeze" their investments to collect the regular coupon payments and benefit from the yields offered by their securities. This strategy is known as “riding the yield curve".

• Assumption 1: The times to maturity of the bonds evaluated at t are assumed to stay constant in $[t; t + \Delta t]$ for $\Delta t$ sufficiently small.

• Assumption 2: The investors are supposed to replicate at t the exact same portfolio constructed at $t-1$ before updating their investments. As soon as $\Delta t$ is such that the earlier approximation does not hold anymore, each aggregated trader sells the corresponding bonds and buys new ones with the appropriate maturities. This replication has no impact on the bonds’ prices as the sell orders of an agent replicating her portfolio are assumed to be compensated by the buy orders of another agent doing the same process. The sole impact on the prices hence results from the excess demands arising from the update of the agents’ investments.

• Assumption 3: The auctions of the coupon paying bonds take place at each time step and are ruled by the Walrasian market making process. They are not reopened, have no impact on the coupon rates and do not distinguish between competitive and non-competitive bidders.

• Assumption 4: The coupon rates are defined by an exogenous stochastic process reflecting the state of the economy. The dependency between the excess demands for the bonds and the coupon rates is ensured by the traders’ investment processes.

• Assumption 5: The outstanding of coupon paying bonds remains constant in time.

Table 1 presents the different types of assets involved in the market model and their corresponding maturities. The preferred stock is granted an infinite maturity and is hence considered as a perpetual bond. The model involves $M = 8$ assets, the first one being the risk-free bond and the last one the preferred stock.

Table 1 Asset types and maturities considered in the market model

The yields to maturity ($\text {YTM}$), implied spot rates, yield duration, convexity and dispersions are defined as follows.

• 1. Yields to maturity The market model is developed such that the bonds’ yields are derived from the associated prices. The prices are updated at each time step through the market making process, with the updates of the corresponding yields being simultaneous. The YTMs are defined as internal rates of returns equalizing the present values of the future cash flows (coupons and principal payments) with the bonds’ prices. They satisfy Eq. (87) presented in Appendix A.1. Their approximate expressions are derived below for each asset $A^k$ with $k \in \{1,\ldots ,M\}$ and where the no-arbitrage prices, times to maturity, principal and coupon payments are, respectively, represented by ${\text {PV}}_k$, $N_k$, $\text {FV}_k$ and $\text {C}_k$:

  $$\begin{aligned} \text {YTM}_k \approx \frac{\text {C}_k + \Big ( \frac{\text {FV}_k - {\text {PV}}_k}{N_k} \Big )}{\frac{\text {FV}_k + {\text {PV}}_k}{2}}. \end{aligned}$$

  (1)

  As the time to maturity of $A^M$ tends to infinity, the associated yield $\text {YTM}_M$ is evaluated by:

  $$\begin{aligned} \text {YTM}_{M} \approx \lim _{N_M\rightarrow \infty } \frac{\text {C}_M + \left( \frac{\text {FV}_M - {\text {PV}}_M}{N_M} \right) }{\frac{\text {FV}_M + {\text {PV}}_M}{2}} = \frac{\text {C}_M}{\frac{\text {FV}_M + {\text {PV}}_M}{2}}. \end{aligned}$$

  (2)

  The consideration of the preferred stock as a perpetual bond is further justified by considering Eq. (89) provided in Appendix A.1, which emphasizes the fact that the price of the associated zero-coupon bond vanishes as the time to maturity tends to infinity.

• 2. Implied spot rates The implied spot rates ($\text {ISR}$) associated with each bond are introduced successively and correspond to the holding period returns realized over the entire life of each bond. They are particular forward rates covering the maturities of the bonds and thus provide approximations of the yields to maturity with the advantage to allow omitting the consideration of the rates at which the coupons are reinvested. These rates are evaluated in terms of the corresponding discount factors, as presented in Appendix A.1.

• 3. Yield durations and convexities The market model considers both modified ($\text {ModDur}^k_t$) and Macaulay ($\text {MacDur}^k_t$) durations associated with each asset. The Macaulay duration is defined as in Marrison (2002). Equations (99) to (101) of Appendix A.1 provide the approximations of these durations as well as of the corresponding convexities ($\text {Conv}^k_t$).

• 4. Portfolio yields, durations, convexities and dispersions The aggregated portfolio statistics are defined successively in Appendix A.1, considering each portfolio as a unique entity equivalent to a bigger bond paying bundles of cash flows at each time step. Two approaches are applied to approximate the yields of the traders’ portfolio: the market-value-weighted (MV-weighted) and basis-point-value-weighted (BPV-weighted) approaches. The latter replaces the traders’ wealth fractions by the proportion of each asset’s BPV to the sum of the BPVs of the assets in which the portfolios are invested. The BPVs are expressed as the change in market value of each asset originated by a variation of its yield by a basis point, i.e. by $0.01 \%$. The modified market-value-weighted average portfolio and Macaulay durations, convexities and dispersions are further defined in the appendices.

This concludes the presentation of the parameters involved in the fixed income framework.

2.2 Assets and wealth dynamics

This section defines the assets’ returns and how the traders’ wealth are updated. A particular attention is paid to the accrued interests in the returns of the coupon-paying bonds of finite maturities.

• Assets’ returns

  • Risk-free asset As in Kaizoji et al. (2015), the risk-free asset provides a constant rate of return $r_f$. $A^1$ is further defined as a zero-coupon bond having a maturity of one year. The selected “day-count” convention defines a year to correspond to 250 time steps. The issuance price ${\text {PV}}^1$ and face value $\text {FV}^1$ of $A^1$ are related by $r_f$ as follows:

    $$\begin{aligned} \text {FV}^1 = {\text {PV}}^1 (1+ 250 \times r_f). \end{aligned}$$

    (3)

  • Preferred stock The preferred stock $A^M$ pays periodical dividends to its holders. Its “one time step” return $r^M_t$ consists of the sum of the price return $\frac{P^M_t}{P^M_{t-1}} -1$ with the return provided by the dividend payment $d^M_t$:

    $$\begin{aligned} \begin{aligned} r^M_t&= \frac{P^M_t + d^M_t}{P^M_{t-1}} -1 \\&= \frac{P^M_t + d^M_{t-1}(1 + r_d + \sigma _d u^M_t)}{P^M_{t-1}}-1, \end{aligned} \end{aligned}$$

    (4)

    where $r_d$ is the long-term growth rate, $\sigma _d$ the standard deviations of the multiplicative process and $u^M_t \sim {\mathcal {N}}(0,1)$ an i.i.d. random variable.

  • Coupon-paying bonds of finite maturities Assuming the absence of taxation and the simultaneous payments of the coupons and dividends, the returns realized on each bond $A^k$ for $k \in \{2,\ldots ,M-1\}$ are expressed as the sum of the price and coupon returns as follows:

    $$\begin{aligned} r^{k}_t = \frac{P^k_t + d^k_t}{P^k_{t-1}} -1. \end{aligned}$$

    (5)

In particular, the coupon payment $d^k_t$ inherently corresponds to the interest accrued between the purchase of the bond $A^k$ at $t-1$ and its sale at t. This arises from the replication process exposed earlier and requiring to buy another bond at t having the same maturity as the bond bought at $t-1$. The accrued interests are related to the annual coupon payments of $C^k_t$ through: $d^k_t = C^k_t \times \mathrm{d}t / 250$, where the time increment is set to $\mathrm{d}t = 1$. The multiplicative process governing the determination of the accrued coupon payments is further expressed as:

$$\begin{aligned} d^k_t = (1 + \sigma ^k \nu ^k_t) d^k_{t-1}, \end{aligned}$$

(6)

where $\sigma ^k$ and $\nu ^k_t$, respectively, represent the standard deviation of the process and the associated i.i.d. random variable $\sim {\mathcal {N}}(0,1)$.

• Wealth dynamics

The wealth dynamics are expressed in terms of the wealth fractions of the corresponding traders and of the assets’ returns described above. Let $W^i_t$ represent the wealth level of trader $i \in \{f,c\}$ at time step t. Also, let the number of shares $z^{i,k}_t$ of asset k held by the corresponding trader at t be expressed in terms of the wealth fraction $x^{i,k}_t$: $z^{i,k}_t = x^{i,k}_t W^{i}_t / P^k_t$. The wealth of the trader allocated to this asset is updated as:

$$\begin{aligned} \Delta W^{i,k}_{t-1 \rightarrow t} = x^{i,k}_{t-1} W^i_{t-1} \frac{P^k_t}{P^k_{t-1}} - x^{i,k}_{t-1} W^i_{t-1}. \end{aligned}$$

(7)

Summing over the 7 coupon- and dividend-paying assets, one obtains:

$$\begin{aligned} W^{i}_{t-1 \rightarrow t} |_{A^2,\ldots ,A^M} = W^i_{t-1} \sum _{k=2}^{M=8}{x^{i,k}_{t-1} r^k_t}. \end{aligned}$$

(8)

The following trivial conditions are moreover applied on the agents’ wealth fractions due to the absence of both borrowing and short selling:

$$\begin{aligned} \left\{ \begin{array}{ll} 0 \le x^{i,k}_t \le 1 \quad \forall t, \\ 0 \le \sum _{k=2}^{M}{x^{i,k}_t} \le 1 \quad \forall t.\end{array}\right. \end{aligned}$$

(9)

The fraction of wealth $x^{i,rf}_t$ invested in the risk-free asset is thus:

$$\begin{aligned} x^{i,rf}_t = 1 - \sum _{k=2}^{M}{x^{i,k}_t}. \end{aligned}$$

(10)

By summing the wealth increments given by Eq. (8) with the one associated with the risk-free asset, one can express the overall wealth dynamics between two consecutive time steps as:

$$\begin{aligned} \Delta W^i_{t-1 \rightarrow t} = W^i_{t-1} \left[ \sum _{k=2}^{M}{x^{i,k}_{t-1} r^k_t} + x^{i,rf}_{t-1} r_f\right] , \end{aligned}$$

(11)

which is equivalent to the formulation proposed hereafter, considering the expressions of the total returns exposed earlier:

$$\begin{aligned} W^i_{t} = W^i_{t-1} \left[ \sum _{k=2}^{M}{x^{i,k}_{t-1} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) } + \left( 1 - \sum _{k=2}^{M}{x^{i,k}_{t-1}} \right) (1+ r_f ) \right] , \end{aligned}$$

(12)

where the coupon and dividend payments associated with each asset $A^k$ for $k \in \{2,\ldots ,8\}$ are:

$$\begin{aligned} d^k_t = \left\{ \begin{array}{ll} (1 + \sigma ^k \nu ^k_t) d^k_{t-1} &{} \text {for}\quad k = 2,\ldots ,7,\\ (1 + r_d + \sigma ^k \nu _t) d^k_{t-1} &{} \text {for}\quad k = 8.\end{array}\right. \end{aligned}$$

(13)

This concludes the presentation of the wealth updates and leads to the introduction of the fundamentalist’s investment strategy.

2.3 Fundamentalist trader

The fundamentalist trader is defined in accordance with the descriptions provided in Chiarella et al. (2009), Kaizoji et al. (2015) and Xu et al. (2014). The trader is granted a constant relative risk-aversion (CRRA) utility function required to be maximized w.r.t. the wealth allocations. This objective is represented by a generalized myopic mean-variance optimization problem expressed in terms of the ex ante expectation of the utility function depending on the asset returns and wealth allocations introduced previously. Accounting for Eq. (12), the optimization problem proposed in Kaizoji et al. (2015) is generalized to the multi-asset case hereafter.

• Generalized optimization problem

$$\begin{aligned} \max _{{\mathbf {x}}^f_t} {\mathrm {E}}_t \left[ U(W^f_{t+1}) \right] \quad \forall t \in [0,T], \end{aligned}$$

(14)

where $U(W) \in {\mathbb {R}}$ is the CRRA utility function and ${\mathbf {x}}^f_t$ the vector containing the wealth fractions associated with each asset at t.

The derivation of the solution is provided in Appendix A.2 and consists in solving the Hamilton–Jacobi–Bellman partial differential equation (PDE) associated with the generalized optimization problem of concern. The equivalent formulation of the optimization problem as a stochastic optimal control problem is proposed below. This is achieved by considering Eq. (12) to define the wealth update in terms of a stochastic differential equation (SDE). To do so, the individual returns $r^k_{t+1}$ of each asset $A^k$ for $k \in \{2,\ldots ,8\}$ are expressed as in Xu et al. (2014) in terms of the expectations ${\mathrm {E}}_t \Big [r^k_{t+1}\Big ]$ formulated at t of the returns to be realized at $t+1$, of the standard deviations $\sigma ^k_t$ and of the i.i.d. random variables $\psi ^k_t$ such that $\psi ^k_t \sim {\mathcal {N}}(0,1)$:

$$\begin{aligned} r^k_{t+1} = {\mathrm {E}}_t \Big [r^k_{t+1} \Big ] + \sigma ^k_{t} \psi ^k_{t}. \end{aligned}$$

(15)

Repeating, we use discrete time equations and the slight abuse of notation of representing one time step by the symbol $\mathrm{d}t$ (which is equal to 1) and the corresponding one-step increment of the wealth of the fundamentalist trader is denoted $\mathrm{d}W^f_t$. Considering Eqs. (12) and (15), the SDE accounting for the stochasticity of the wealth level of the fundamentalist is hence expressed as:

$$\begin{aligned} \begin{aligned} \mathrm{d}W^f_t =&W^f_{t-1} \left[ \left( 1 - \sum _{k=1}^{M-1}{x^{f,k}_{t-1}} \right) r_f + \sum _{k=1}^{M-1}{x^{f,k}_{t-1}}{\mathrm {E}}_{t-1} \left[ r^k_{t}\right] \right] \\&\quad + W^f_{t-1} \sum _{k=1}^{M-1}{x^{f,k}_{t-1} \sigma ^k_{t-1} \psi ^k_{t-1}}, \end{aligned} \end{aligned}$$

(16)

which can be synthesized by the successive formulation proposed in Xu et al. (2014),

$$\begin{aligned} \mathrm{d}W^f_t = \mu _{t}(W^f) \mathrm{d}t + \sigma _{t} (W^f) \Psi _t, \end{aligned}$$

(17)

where $\mu _t$ represents the drift coefficient, $\sigma _t$ the standard deviation and $\Psi _t \sim {\mathcal {N}}(0,1)$ an i.i.d. random variable. Considering the unit vector ${\mathbf {e}} \in {\mathbb {R}}^{M-1}$ and the transpose operator $(.)^\top$, one can further express $\mu _t (W^f)$ and $\sigma _t (W)$ in both scalar and vector notations:

$$\begin{aligned} \begin{array}{ll} \mu _t (W^f) &{}= W^f_{t-1} \left[ \left( 1 - \sum _{k=1}^{M-1}{x^{f,k}_{t-1}} \right) r_f + \sum _{k=1}^{M-1}{x^{f,k}_{t-1}} {\mathrm {E}}_{t-1} \left[ r^k_{t}\right] \right] \\ &{} = W^f_{t-1} \left[ \left( 1-{{\mathbf {x}}^f}_{t-1}^\top {\mathbf {e}}\right) r_f + {{\mathbf {x}}^f}_{t-1}^\top {\mathrm {E}}_{t-1} \left[ {\mathbf {r}}_{t}\right] \right] , \end{array} \end{aligned}$$

(18)

$$\begin{aligned} \begin{array}{ll} \sigma _t(W^f) &{} = \left\{ {W^f_{t-1}}^2 \left[ \sum _{k=1}^{M-1}{{x^{f,k}_{t-1}}^2 {\sigma ^{k}_{t-1}}^2 + \sum _{k=1}^{M-1}{\sum _{l=1}^{M-1}{x^{f,k}_{t-1}x^{f,l}_{t-1} \sigma ^{kl}_{t-1}}}}\right] \right\} ^\frac{1}{2} \\ &{} = W^f_{t-1} {{\mathbf {x}}^f}_{t-1}^\top {\sigma }_{t-1}, \end{array} \end{aligned}$$

(19)

where we distinguish $\sigma _t (W^f)$ and ${\sigma }_t$, the former accounting for the diffusion term of Eq. (17) and the latter representing the vector containing the standard deviations of the individual assets’ returns. Combining Eqs. (17), (18) and (19), the resulting SDE is given by:

$$\begin{aligned} \mathrm{d}W^f_t = W^f_{t-1} \left[ \left( 1-{{\mathbf {x}}^f}_{t-1}^\top {\mathbf {e}}\right) r_f + {{\mathbf {x}}^f}_{t-1}^\top {\mathrm {E}}_{t-1} \left[ {\mathbf {r}}_{t}\right] \right] \mathrm{d}t + W^f_{t-1} {{\mathbf {x}}^f}_{t-1}^\top {\sigma }_{t-1} \Psi _t, \end{aligned}$$

(20)

where $\Psi _t$ represents a one-dimensional Brownian motion and where the covariance matrix is expressed as $\Sigma = \sigma \sigma ^T \in {\mathbb {R}}^{(M-1) \times (M-1)}$. Letting $\mu \in {\mathbb {R}}^{M-1}$ account for ${\mathrm {E}}_{t-1} \left[ {\mathbf {r}}_{t}\right]$, one can finally write:

$$\begin{aligned} \mathrm{d}W^f = W^f \Big [(1-{\mathbf {x}}^\top {\mathbf {e}})r_f + {\mathbf {x}}^\top \mu \Big ] \mathrm{d}t + W^f {\mathbf {x}}^\top {\sigma } \Psi . \end{aligned}$$

(21)

The cost functional subjected to the previous wealth dynamics is subsequently introduced as follows:

$$\begin{aligned} {\mathcal {J}}({\mathbf {x}}) = {\mathrm {E}}[U(W^f)] \quad \in {\mathbb {R}}. \end{aligned}$$

(22)

The problem hence requires to find the optimal control ${\mathbf {x}} \in {\mathbb {R}}^{M-1}$ maximizing the cost functional subjected to the wealth defined by the stochastic system of Eq. (21).

• Generalized stochastic optimal control problem:

$$\begin{aligned} {\mathcal {J}}({\mathbf {x}}^*) = \max _{{\mathbf {x}}:[0,T] \mapsto {\mathbb {R}}^{M-1}} {\mathcal {J}}({\mathbf {x}}). \end{aligned}$$

(23)

The solution is given by the optimal state trajectory ${\mathbf {x}}^* : [0,T] \mapsto {\mathbb {R}}^{M-1}$. Introducing the risk aversion parameter $\gamma \ne 1$ ($\gamma > 0$ always holds), the cost-to-go function is given by:

$$\begin{aligned} {\mathcal {J}}({\mathbf {x}},t) = {\mathrm {E}}\left[ \frac{W^f(t)^{1-\gamma }}{1-\gamma }\right] . \end{aligned}$$

(24)

Consequently, the stochastic optimal control problem is expressed as follows, considering the differential form of the standard Brownian motion $\Theta (t) \in {\mathbb {R}}$:

$$\begin{aligned} \max _{{\mathbf {x}} : [0,T] \mapsto {\mathbb {R}}^{M-1}} {\mathrm {E}}\left[ \frac{W^f(t)^{1-\gamma }}{1-\gamma }\right] \end{aligned}$$

(25)

subjected to:

$$\begin{aligned} \begin{aligned} \mathrm{d}W^f =&\mu (W^f)\mathrm{d}t + \sigma (W^f)\mathrm{d}\Theta \\ W^f(0) =&W^f_0. \end{aligned} \end{aligned}$$

(26)

The solution of the generalized stochastic optimal control problem is obtained in Appendix A.2 and reads

$$\begin{aligned} {x_k}^* = \frac{1}{\gamma } \sum _{l=1}^{M-1}{ {\Sigma }^{-1}_{kl} } (\mu _{l} - r_f) \quad \text {for} \quad k \in \{2,\ldots ,M\}. \end{aligned}$$

(27)

The remaining unknowns consist of the covariance matrix and of the vector containing the expectations of the future returns. The matrix ${{\varvec{\Sigma }}} = {\sigma } {\sigma }^\top$ is defined to be diagonal and the individual variances remain constant. One obtains the following notations introducing the Kronecker delta $\delta _{kl}$:

$$\begin{aligned} \left\{ \begin{array}{ll} {\Sigma }_{kl} = \delta _{kl} {\sigma _{kl}}^2, \\ {\Sigma }_{kl}^{-1} = \frac{\delta _{kl}}{{\sigma _{kl}}^2}. \end{array}\right. \end{aligned}$$

(28)

The expectations of the future returns are taken as the sum of a constant and a variable term for each of the yields to maturity and implied spot rates. This allows us to account for the trader’s own opinion about the returns provided by each bond (constant term) while still considering the reality of the market (variable term). The constant terms consist of the initial values of the parameters. The variable terms are those evaluated at $t-1$ as a consequence of the replication process explained in Sect. 2.1. Therefore, for each bond $A^k$ where $k \in \{2,\ldots ,M-1\}$, one expresses:

$$\begin{aligned} \mu ^k_t= & {} \omega _1 \frac{1}{2} \Big \{ \Big [(1 +\text {YTM}^k_0)^{\frac{1}{250}}-1 \Big ] + \Big [(1 + \text {ISR}^k_0)^{\frac{1}{250}}-1 \Big ] \Big \} + \omega _2 \Big [ (1+\text {YTM}^k_{t-1})^{\frac{1}{250}} -1\Big ] \nonumber \\&+ \omega _3 \Big [ (1+ \text {ISR}^k_{t-1})^{\frac{1}{250}} -1 \Big ] + \frac{d^k_t}{P^k_t}, \end{aligned}$$

(29)

where $\omega _1$, $\omega _2$ and $\omega _3$ are the weights, respectively, associated with the initial values of the yields to maturity and implied spot rates as well as with their evaluations at $t-1$. As explained in Appendix A.1, the implied spot rates cannot be computed for the preferred stock. Moreover, due to the particular identity of the asset, the varying term associated with the yield to maturity is neglected. For $k = M$, one hence obtains:

$$\begin{aligned} \mu ^M_t = (1 +\text {YTM}^M_0)^{\frac{1}{250}}-1 + \frac{d^M_t(1 + r_d)}{P^M_t}. \end{aligned}$$

(30)

The wealth fractions of the fundamentalist trader are ultimately expressed as:

$$\begin{aligned} x^{f,k}_t = \frac{1}{\gamma \sigma ^2_k} (\mu ^k_t - r_f) \quad \text {for} \quad k \in \{2,..,M\}, \end{aligned}$$

(31)

with the boundary condition given by $\sum _{i=1}^M x^{f,i}_t=1$.

The subsequent developments introduce the investment process of the chartist trader taking the other side of the trades.

2.4 Noise traders

As in Kaizoji et al. (2015), the individual chartist traders do not diversify their allocations but instead select an individual asset in which they invest their entire wealth. Their investments are updated probabilistically at each time step through an Ising model-like set-up under the influence of the momentum in the assets’ prices’ time series and the opinion of other noise traders. At the aggregated level, the representative noise trader allocates her wealth proportionally to the number of individual traders invested in the corresponding assets.

Fig. 1

Representation of the market as seen by the chartists: differentiation of the two pools of assets and corresponding steps of the investment process. The arrows represent the transit of the agents’ wealth

As shown in Fig. 1, the noise trader identifies two mutually exclusive pools of assets. The first one contains all the bonds including the risk-free one and is identified as ${\mathcal {P}}$. The second one is solely constituted of the preferred stock $A^M$. Built on this market representation, three distinct steps are defined below to decompose the investment process of the aggregated chartist. Figure 1 further illustrates which particular pool is concerned by each of these steps. The first one tackles the rebalancing of the aggregated chartist’s portfolio between the two asset classes considered. The second one is concerned with the update of the investments of the individual agents inside the bonds of ${\mathcal {P}}$. The third one finally reconciles the two previous updates simultaneously in order to derive the wealth fractions of the aggregated trader.

Step 1: Aggregated chartist’s portfolio rebalancing between ${\mathcal {P}}$ and $A^M$:

The first step of the investment process of the aggregated chartist deals with the update of the allocations in either of the two groups of assets defined previously. It derives an estimate of the number of individual traders transiting between the two asset classes. Both social imitation and adherence to a momentum-following strategy influence the outcome of the present step.

The number of individual chartists invested in $A^{M}$ is $N^{M}$ and the number of remaining agents is $N^{1 \circlearrowleft M-1}$. The total fraction of wealth invested in $A^M$ at t is hence:

$$\begin{aligned} \left. x^c_t \right| _{A^{M}} = \frac{N^{M}_t}{N^{M}_t+N^{1 \circlearrowleft M-1}_t}. \end{aligned}$$

(32)

Following the definition of Kaizoji et al. (2015), the associated opinion index $s_t^c$ is defined as:

$$\begin{aligned} s^{\text {step}1}_t = \frac{N^{M}_t - N^{1 \circlearrowleft M-1}_t}{N^{M}_t + N^{1 \circlearrowleft M-1}_t} = 2 \left. x^c_t \right| _{A^{M}} -1 \quad \in [-1,1]. \end{aligned}$$

(33)

The parameter accounting for momentum is moreover defined as:

$$\begin{aligned} H^{\text {step1}}_t = H^M_t - H^{1 \circlearrowleft M-1}, \end{aligned}$$

(34)

introducing the momentum associated with the preferred stock:

$$\begin{aligned} H^M_t = \theta H^M_{t-1} + (1 - \theta ) \left( \frac{P^M_{t}}{P^M_{t-1}} -1 \right) \end{aligned}$$

(35)

and the one related to ${\mathcal {P}}$ as:

$$\begin{aligned} H^{1 \circlearrowleft M-1}_t = \theta H^{1 \circlearrowleft M-1}_{t-1} + (1 - \theta ) \left[ (1+\text {YTM}^{M-1}_{t})^{\frac{1}{250}} - (1+\text {YTM}^{1}_{t})^{\frac{1}{250}} \right] . \end{aligned}$$

(36)

Parameter $\theta$ controls the time scale $\sim {1 \over 1-\theta }$ over which the momentum is estimated.

Equation (36) introduces the simple “rule of thumb” followed by the individual chartists to compare the price returns of the isolated asset with the “term-spreads” provided by ${\mathcal {P}}$. As shown in Eq. (34), this spread is defined as the difference between the yields to maturity of the bonds situated at the extremities of the yield curve, fitted to the “one-time-step” returns of the preferred stock. This spread moreover represents a shorthand evaluation of the steepness of the yield curve giving insights about potential future returns obtained in excess of those given by the risk-free bond. Considering the fact that $\text {YTM}^{1}_{t} = \text {cst.}$, the higher the yields of the longer-term bond, the more the aggregated chartist is influenced to allocate her wealth in the fixed income asset class.

As exposed in Kaizoji et al. (2015), the herding propensity $\kappa _t$ is introduced in order to represent the noise traders’ susceptibility to herding and their propensity to adhere to a momentum following strategy. This parameter is moreover defined to be either constant or to follow a stochastic mean-reverting process. In the latter case, it incorporates the alternating regimes of pessimistic mood and exuberance observed in financial markets, as detailed in Shiller (2015) and in Sornette (2003). In the former case, it steers the model dynamics by accounting for stable market regimes. Section 2.6 provides a thorough definition of the two distinct processes ruling the dynamics of the parameter. Given the aforementioned rationale for its introduction in the model, the herding propensity is directly impacting the noise traders’ dynamics through the transition probabilities presented in Eq. 37 and involving the opinion index and momentum parameter.

By definition, the transition probabilities $p^+_t$ and $p^-_t$ are associated with the following actions: an individual noise trader holding the preferred stock $A^M$ at t shifts her investment to ${\mathcal {P}}$ with a probability $p_t^+$ and another noise trader invested in ${\mathcal {P}}$ at the same moment transits to $A^M$ with a probability $p_t^-$. These probabilities are consequently expressed as:

$$\begin{aligned} p^{\pm }_t(s^{\text {step}1}_t,H^{\text {step1}}_t) = \frac{1}{2}p^{\pm } \left[ 1 \mp \kappa _t \left( s^{\text {step}1}_t + H^{\text {step1}}_t \right) \right] , \end{aligned}$$

(37)

where $p^+ = \text {cst.}$ and $p^- = \text {cst.}$ control the average holding time associated with each group of assets. Setting $p^- > p^+$ implies that an up- or downward increment of either of $H^{\text {step1}}_t$ or $s^{\text {step}1}_t$ is not engendering the same reaction of the agents. This is translated by a favour given to $A^M$ to the detriment of ${\mathcal {P}}$.

Introducing the Bernoulli random variables $\xi _k(p)$ taking the value of 1 with a probability of p and 0 otherwise, one can express the subsequent master equations:

$$\begin{aligned} \left\{ \begin{array}{ll} N^{M}_t = \sum _{k=1}^{N^{M}_{t-1}}{ [1 - \xi _k(p^+_{t-1})]} + \sum _{k=1}^{N^{1 \circlearrowleft M-1}_{t-1}}{ \xi _k(p^-_{t-1})}, \\ N^{{1 \circlearrowleft M-1}}_t = \sum _{k=1}^{N^{1 \circlearrowleft M-1}_{t-1}}{ [1 - \xi _k(p^-_{t-1})]} + \sum _{k=1}^{N^{M}_{t-1}}{ \xi _k(p^+_{t-1})}. \end{array}\right. \end{aligned}$$

(38)

The net number of agents transiting between both pools of assets at t is finally estimated by the absolute value of the variable $\Delta N^{+ \rightarrow -}_t$ defined as:

$$\begin{aligned} \left|\Delta N^{+ \rightarrow -}_t \right|= \left|\sum _{k=1}^{N^{M}_{t-1}}{ \xi _k(p^+_{t-1})} - \sum _{k=1}^{N^{1 \circlearrowleft M-1}_{t-1}}{ \xi _k(p^-_{t-1})} \right|. \end{aligned}$$

(39)

Step 2: Aggregated chartist’s investments within ${\mathcal {P}}$:

The second step of the allocation process tackles the allocations of the chartists inside ${\mathcal {P}}$. Their updates are achieved independently for each pair of neighbouring assets $(A^i, A^{i+1})$ and $(A^1, A^{M-1})$. Again, both social imitation and adherence to a momentum-following strategy influence the outcome of these updates. The opinion index associated with the pair $(A^i, A^{i+1})$ is defined below, fitting the variable to the momenta introduced subsequently and introducing $N^i_t$ as the number of chartists invested in $A^i$ at t:

$$\begin{aligned} \left. s^{\text {step2}}_t \right| _{i,i+1} = \frac{1}{2} \frac{N^{i+1}_t - N^i_t}{N^{i+1}_t + N^i_t}. \end{aligned}$$

(40)

The momentum parameter associated with the same pair of assets is moreover expressed as:

$$\begin{aligned} \begin{aligned} \left. H^{\text {step2}}_t \right| _{i,i+1}&= H^{i+1}_t - H^i_t, \end{aligned} \end{aligned}$$

(41)

introducing the momentum parameters associated with each asset. These are further expressed in terms of the sub-parameters accounting for the momenta observed in the yields and in the prices:

$$\begin{aligned} H^i_t = \alpha H^{i,\text {price}}_t + (1-\alpha ) H^{i,\text {yield}}_t, \end{aligned}$$

(42)

where $\alpha \in [0,1]$. In particular:

$$\begin{aligned} H^{i,\text {price}}_t= & {} \left\{ \begin{array}{ll} \theta H^{i,\text {price}}_{t-1} + (1- \theta ) \left( \frac{P^i_t}{P^i_{\text {auction}}} -1 \right) &{} \text {for}\quad i = 2,\ldots ,M-1, \\ r_f&{} \text {for}\quad i = 1. \end{array}\right. \end{aligned}$$

(43)

$$\begin{aligned} H^{i,\text {yield}}_t= & {} \theta H^{i,\text {yield}}_{t-1} + (1-\theta ) \Big [ (1+\text {YTM}^i_{t-1})^{\frac{1}{250}} - (1+\text {YTM}^i_{t})^{\frac{1}{250}} \Big ], \end{aligned}$$

(44)

where the yields are fitted to the “one time step” periods over which the price returns are determined. Equation (44) is moreover introducing the consideration of decreasing yields by the chartists. This enforces the assessment of positive realized price returns and contrasts with the previous step in which the agents valued higher yields implying higher potential price returns when looking at the overall market. Once invested inside the fixed income asset class, these agents naturally favour decreasing yields corresponding to positive realized price returns.

In accordance with the formulation of Eq. (37), the transition probabilities associated with each pair of consecutive assets are defined as:

$$\begin{aligned} {p^{i,i+1}_t}^{\pm } = \frac{p^{\pm }}{2} \Big [1 \mp \kappa _t \Big (\left. s^{\text {step2}}_t \right| _{i,i+1} + \left. H^{\text {step2}}_t \right| _{i,i+1}\Big )\Big ], \end{aligned}$$

(45)

where ${p^{i,i+1}_t}^{+}$ quantifies the probability to invest in $A^{i}$ when starting from $A^{i+1}$. Accordingly, ${p^{i,i+1}_t}^{-}$ accounts for the probability to invest in $A^{i+1}$ when starting from $A^{i}$.

The number of agents transiting between two neighbouring assets in either direction is hence:

$$\begin{aligned} \left\{ \begin{array}{ll} \Delta N^{{i+1 \rightarrow i}}_t = \sum _{k=1}^{N^{i+1}_{t-1}}{ \xi _k \left( {p^{i,i+1}_{t-1}}^+\right) }, \\ \Delta N^{{i \rightarrow i+1}}_t = \sum _{k=1}^{N^{i}_{t-1}}{ \xi _k \left( {p^{i,i+1}_{t-1}}^-\right) }. \end{array}\right. \end{aligned}$$

(46)

Ultimately, the net number of agents leaving or entering each bond i is expressed as:

$$\begin{aligned} \left. \Delta N^{\text {step2}}_t \right| _{i} = \Delta N^{{i-1 \rightarrow i}}_t - \Delta N^{{i \rightarrow i-1}}_t + \Delta N^{{i+1 \rightarrow i}}_t -\Delta N^{{i \rightarrow i+1}}_t . \end{aligned}$$

(47)

Step 3: Execution of the update of the aggregated chartist’s investments:

The last step executes the update of the investments on the aggregated level. The following opinion index is first introduced for each of the $M-1$ bonds of fixed maturities:

$$\begin{aligned} \left. s^{\text {step3}}_t \right| _{i} = \frac{(M-1) N^i_t - \sum _{k=1, k \ne i}^{M-1}{N^k_t}}{(M-1)N^i_t + \sum _{k=1, k \ne i}^{M-1}{N^k_t}}. \end{aligned}$$

(48)

The performance indicators accounting for the popularity of each bond and for the trend observed in their returns are defined as follows, considering the parameters $H^i_t$ introduced in Eq. (42):

$$\begin{aligned} \text {Perf}^i_t = \left. s^{\text {step3}}_t \right| _{i} + H^i_t. \end{aligned}$$

(49)

Building on these parameters, the categorical probabilities are further defined as

$$\begin{aligned} \left. p^{\text {step3}}_t \right| _{i} = \left\{ \begin{array}{ll}\frac{\text {Perf}^i_t - \text {Perf}^{\text {worst}}_t}{\sum _{k=1}^{M-1}{( \text {Perf}^k_t - \text {Perf}^{\text {worst}}_t )}} &{} \text {if}\quad \Delta N^{+ \rightarrow -}_t > 0, \\ \frac{\text {Perf}^{\text {best}}_t - \text {Perf}^i_t}{\sum _{k=1}^{M-1}{(\text {Perf}^{\text {best}}_t - \text {Perf}^k_t)}} &{} \text {if}\quad \Delta N^{+ \rightarrow -}_t < 0,\end{array}\right. \end{aligned}$$

(50)

where $\text {Perf}^{\text {worst}}_t$ and $\text {Perf}^{\text {best}}_t$ are, respectively, $\min _{i \in \{1,\ldots ,M-1\}} \text {Perf}^{i}_t$ and $\max _{i \in \{1,\ldots ,M-1\}} \text {Perf}^{i}_t$.

The weights $\omega ^i_t$ associated with each bond are finally defined as:

$$\begin{aligned} \omega ^i_t = \frac{N^i_t \left. p^{\text {step3}}_t \right| _{i}}{\sum _{k=1}^{M-1}{ N^k_t \left. p^{\text {step3}}_t \right| _{k}}}. \end{aligned}$$

(51)

As a matter of fact, $\omega ^i_t$ corresponds to the share of the expected value of the Bernoulli trial associated with $A^i$ and endowed with a probability of success of $\left. p^{\text {step3}}_t \right| _{i}$ over the sum of each of the $M-1$ expectations. These weights are introduced to adjust the dispersion of the agents inside ${\mathcal {P}}$ as can be inferred from the expression of Eq. (52) presented below.

The number of agents transiting either to or from each bond of ${\mathcal {P}}$ is subsequently introduced as $\left. \Delta N^{\text {step3}}_t \right| _{i}$. By considering the result of the update achieved in the first step of the allocation process, one can hence write:

$$\begin{aligned} \left. \Delta N^{\text {step3}}_t \right| _{i} = \big \lfloor \omega ^i_t \Delta N^{+ \rightarrow -}_t \big \rfloor , \end{aligned}$$

(52)

where $\lfloor . \rfloor$ represents the floor operator dropping the decimal part. With this variable, one can execute the update of the individual investments considering the outcome of the two previous steps:

$$\begin{aligned} N^i_t = N^i_{t-1} + \left. \Delta N^{\text {step2}}_t \right| _{i} + \left. \Delta N^{\text {step3}}_t \right| _{i} \quad \ge 0 \quad \forall t. \end{aligned}$$

(53)

Accordingly, the updated number of agents invested in the preferred stock is given by:

$$\begin{aligned} N^M_t = N^M_{t-1} - \sum _{k=1}^{M-1}{\left. \Delta N^{\text {step3}}_t \right| _{i} }. \end{aligned}$$

(54)

As a result, the updated share of wealth of the aggregated chartist trader invested in each asset $A^i$ for $i \in \{1,\ldots ,M\}$ is ultimately expressed as:

$$\begin{aligned} \left. x^c_t \right| _{i} = \frac{N^k_t}{\sum _{k=1}^{M}{N^k_t}}. \end{aligned}$$

(55)

This concludes the presentation of the aggregated noise trader and leads to the formulation of the market making process.

2.5 Market clearing process

The present section tackles the market clearing process executed at each time step by the Walrasian auctioneer. First, the expressions of the excess demands formulated by both aggregated traders are developed. The system of nonlinear equations resulting from the satisfaction of the Walrasian equilibrium condition is provided subsequently.

On the one hand, the excess demand of the fundamentalist trader $\Delta D^{f,i}_{t-1 \rightarrow t}$ associated with each asset $A^i$ for $i \in \{2,\ldots ,M\}$ is expressed at each t as:

$$\begin{aligned} \Delta D^{f,i}_{t-1 \rightarrow t}= & {} W^f_{t-1} \left( A^i_t + \frac{B^i_t}{P^i_t} \right) \left[ \left( 1 - \sum _{k=2}^{M}{x^{f,k}_{t-1}} \right) (1+r_f) + \sum _{k=2}^{M}{x^{f,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] \nonumber \\&- x^{f,i}_{t-1} W^f_{t-1} \frac{P^i_t}{P^i_{t-1}}, \end{aligned}$$

(56)

where $A^i_t$ and $B^i_t$ are obtained from Eq. (31) presenting the expression of the wealth fractions associated with each asset:

$$\begin{aligned} A^i_t= & {} \frac{\Omega ^i_t - r_f}{\gamma \Sigma _{ii}}, \end{aligned}$$

(57)

$$\begin{aligned} B^i_t= & {} \left\{ \begin{array}{ll} \frac{d^i_t}{\gamma \Sigma _{ii}} &{} \text {for}\quad i = 2,\ldots ,M-1, \\ \frac{d^M_t (1 + r_d)}{\gamma \Sigma _{MM}} &{} \text {for}\quad i = M.\end{array}\right. \end{aligned}$$

(58)

As a reminder and considering Eq. (29), $\Omega ^i_t$ for $i \in \{2,\ldots ,M-1\}$ is further expressed as:

$$\begin{aligned} \Omega ^i_t= & {} \omega _1 \frac{1}{2} \left\{ \left[ (1 +\text {YTM}^i_0)^{\frac{1}{250}}-1 \right] + \left[ (1 + \text {ISR}^i_0)^{\frac{1}{250}}-1 \right] \right\} \nonumber \\&+ \omega _2 \left[ (1+\text {YTM}^i_{t-1})^{\frac{1}{250}} -1\right] + \omega _3 \left[ (1+ \text {ISR}^i_{t-1})^{\frac{1}{250}} -1 \right] . \end{aligned}$$

(59)

And for $i = M$, $\Omega ^M_t$ is:

$$\begin{aligned} \Omega ^M_t = (1 +{\text {YTM}}^M_0)^{\frac{1}{250}}-1. \end{aligned}$$

(60)

On the other hand, the excess demand of the aggregated chartist trader associated with the same asset is expressed as:

$$\begin{aligned} \Delta D^{c,i}_{t-1 \rightarrow t}= & {} W^c_{t-1} x^{c,i}_{t} \left[ \left( 1 - \sum _{k=2}^{M}{x^{c,k}_{t-1}} \right) (1+r_f) + \sum _{k=2}^{M}{x^{c,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] \nonumber \\&- x^{c,i}_{t-1} W^c_{t-1} \frac{P^i_t}{P^i_{t-1}}. \end{aligned}$$

(61)

As has been exposed in the second step of the allocation process of the chartist trader presented in Sect. 2.4, the returns of the individualized bonds are evaluated from the auction prices. These being further set to $P^i_{\text {auction}} =1$, the return of each bond realized at t is hence $P^i_t / P^i_{\text {auction}} - 1 = P^i_t - 1$. $A^M$ being further exempted from any auction, its returns are directly evaluated by $P^M_t / P^M_{t-1}$. As a consequence, the following rule is applied in Eq. (61):

$$\begin{aligned} P^i_{t-1} = \left\{ \begin{array}{ll} 1 &{} \text {for}\quad i \in \{2,\ldots ,M-2\}, \\ P^i_{t-1} &{} \text {for}\quad i = M.\end{array}\right. \end{aligned}$$

(62)

For all $i \in \{2,\ldots ,M\}$, the Walrasian equilibrium condition is finally expressed as:

$$\begin{aligned} \Delta D^{f,i}_{t-1 \rightarrow t} + \Delta D^{c,i}_{t-1 \rightarrow t} = 0. \end{aligned}$$

(63)

The system of nonlinear equations involving the unknown prices remaining to be evaluated at each time step is finally obtained by inserting Eqs. (56) and (61) in Eq. (63). The resulting system is presented below and the detailed developments leading to this formulation are provided in Appendix A.3. The equation determining the prices reads

$$\begin{aligned} {P^i_t}^2 (\alpha _{ii} - \beta _i ) + P^i_t (\zeta _i + \chi _{ii}) + \sum _{k=2, k \ne i}^{M}{P^k_t P^i_t \alpha _{ik}} + \sum _{k=2, k \ne i}^{M}{P^k_t \chi _{ik}} + \lambda _i = 0, \end{aligned}$$

(64)

where ${\alpha } \in {\mathbb {R}}^{(M-1) \times (M-1)}$, ${\beta } \in {\mathbb {R}}^{M-1}$, ${\zeta } \in {\mathbb {R}}^{M-1}$, ${\chi } \in {\mathbb {R}}^{(M-1) \times (M-1)}$ and ${\lambda } \in {\mathbb {R}}^{M-1}$ are defined as:

$$\begin{aligned} \left\{ \begin{array}{lllll} \alpha _{ij} = \frac{x^{f,j}_{t-1} W^f_{t-1} A^i_{t} + x^{c,j}_{t-1} W^c_{t-1} x^{c,i}_{t-1}}{P^j_{t-1}}, \\ \beta _i = \frac{x^{f,i}_{t-1} W^f_{t-1} + x^{c,i}_{t-1} W^c_{t-1}}{P^i_{t-1}}, \\ \zeta _i = W^f_{t-1} A^i_{t} \left( 1 - \sum _{k=2}^{M}{x^{f,k}_{t-1}}\right) (1+r_f) + W^c_{t-1} x^{c,i}_t \left( 1 - \sum _{k=2}^{M}{x^{c,k}_{t-1}}\right) (1+r_f) \\ \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad \quad + \sum _{k=2}^{M}{\left[ \frac{d^k_t}{P^k_{t-1}} \left( x^{f,k}_{t-1} W^f_{t-1} A^i_{t} + x^{c,k}_{t-1} W^c_{t-1} x^{c,i}_{t} \right) \right] }, \\ \chi _{ij} = \frac{x^{f,j}_{t-1} B^i_{t}}{P^j_{t-1}}, \\ \lambda _i = W^f_{t-1} B^i_{t} \left[ \left( 1 - \sum _{k=2}^{M}{x^{f,k}_{t-1}}\right) (1+r_f) + \sum _{k=2}^{M}{x^{f,k}_{t-1} \frac{d^k_t}{P^k_{t-1}}} \right] . \end{array}\right. \end{aligned}$$

(65)

Ultimately, this system of $M-1$ equations can be synthesized hereafter by the functions $f^i(t) \in {\mathbb {R}}$:

$$\begin{aligned} f^i (P^2_t,..,P^i_t,\ldots ,P^{M-1}_t) = 0 \quad \text {for} \quad i \in \{2,M-1\}. \end{aligned}$$

(66)

An initial guess of the solution is provided by the previously realized prices $(P^2_{t-1},\ldots ,P^i_{t-1},\ldots P^{M-1}_{t-1})$ and is further refined by the following process introducing ${\mathbf {P}}$, the vector notation of the prices of the coupon and dividend paying assets:

$$\begin{aligned} {\mathbf {P}} \rightarrow {\mathbf {P}}' = {\mathbf {P}} - {\varvec{{{\mathcal {J}}}}}^{-1}{\mathbf {f}}({\mathbf {P}}), \end{aligned}$$

(67)

where ${\varvec{{{\mathcal {J}}}}}$ is the Jacobian matrix defined by ${\mathcal {J}}_{ij} = \partial f^i / \partial P^j$. This concludes the presentation of the market clearing process.

2.6 Dynamics of the herding propensity $\kappa _t$

Recall that the herding propensity $\kappa _t$ represents the noise traders’ susceptibility to herding and their propensity to adhere to a momentum following strategy. It is defined such as to have a direct impact on the transition probabilities $p^+_t$ and $p^-_t$ given by expression (37) in Sect. 2.4. Following Kaizoji et al. (2015), we endow $\kappa _t$ with a dynamics capturing the changing nature of economic and social conditions that influence the behaviour of noise traders. As economic and social conditions are arguably not systematically constant, it is reasonable to consider market regimes corresponding to $\kappa \ne cst.$, thereby introducing a source of stochastic variability inside the transition probabilities through the varying herding propensity. On the other hand, economic and social conditions happen to evolve relatively slowly during certain periods and are characterized by tendencies, trends and variability that can be considered approximately stationary. This scenario is accounted for by $\kappa = cst.$, where the influence of the herding propensity on the transition probabilities remains steady. In fact, fixing the herding propensity to a constant value is analogous to neutralizing the exogenous influence of the social and economic conditions, thereby rendering the model a closed autonomous system. Hence, we will mainly consider the case where $\kappa _t$ follows a stochastic mean-reverting process defined as

$$\begin{aligned} \kappa _t = \kappa _{t-1} + \eta _k (\mu _k - \kappa _{t-1}) + \sigma _k \nu _t, \end{aligned}$$

(68)

where the initial condition is given by $\kappa _0$ and where $\nu _t$ is an i.i.d. discrete-time white noise process with a mean of 0 and a standard deviation of 1. The mean-reversion rate is further represented by $\eta _k$, the mean value by $\mu _k$ and the diffusion associated with the Wiener process $\nu _t$ by $\sigma _k$.

Table 2 Expectation and variance of the herding propensity $\kappa _t$

The instantaneous and long-run evaluations of the expectation and variance of the herding propensity are successively given in Table 2 above. The process of concern is hence found to be stationary and normally distributed in the long-run, such that the following expression holds:

$$\begin{aligned} \kappa _t \sim {\mathcal {N}} \left( \mu _k,\frac{\sigma _k^2}{2 \eta _k} \right) . \end{aligned}$$

(69)

The statistical estimation of the parameters of Eq. 68 is finally achieved as in Kaizoji et al. (2015) via the following expressions:

$$\begin{aligned} \eta _k= & {} \frac{1}{\Delta T} \log \left( \frac{0.2 }{1 - \mu _k} \right) , \end{aligned}$$

(70)

$$\begin{aligned} \sigma _k= & {} 0.1 \sqrt{2 \eta _k}, \end{aligned}$$

(71)

where $\Delta T$ is the time frame during which $\kappa _t$ is set to revert when being in the super-critical regime identified by the bounds located two standard deviations away from the mean.

3 Model dynamics

3.1 Parameter selection

The initialization of the parameters introduced in the market model is subsequently exposed hereafter. The analyses of the time series generated in typical simulations ensue.

The wealth fractions of the aggregated traders are first initialized below, starting with the fundamentalist trader. Given the initial wealth fraction $x^{f,M}_0$, the constant variance $\sigma _M^2$ and the initial expectation of the future return $\mu ^M_0 = {\mathrm {E}}^M_{r,t} + d_0 (1+r_d) / P_0$ associated with $A^M$, the constant risk aversion parameter $\gamma$ is defined a priori by the following expression:

$$\begin{aligned} \gamma = \frac{\mu ^M_0 - r_f}{x^{f,M}_0 \sigma ^2_M}. \end{aligned}$$

(72)

The expressions of the initial wealth fractions associated with the coupon paying bonds $A^i$ for $i \in \{2,\ldots ,M-1\}$ are obtained accordingly:

$$\begin{aligned} x^{f,i}_0 = \frac{\mu ^i_0 - r_f}{\gamma \sigma ^2_i}, \end{aligned}$$

(73)

by ensuring that the following condition is satisfied:

$$\begin{aligned} 0 \le \sum _{k=2}^{M}{x^{f,i}_t} \le 1 \quad \text {for}\quad t = 0. \end{aligned}$$

(74)

The expectations of the future returns $\mu ^i_0$ associated with these assets are defined below and the constant values of the variances $\sigma ^2_i$ are provided in Table 3.

On the other hand, the initial wealth fractions of the aggregated chartist are obtained from the value of $x^{c,M}_0$ through the following expression:

$$\begin{aligned} x^{c,i}_0 |_{i \ne 1,M} = \frac{1 - x^{c,M}_0}{M-1}. \end{aligned}$$

(75)

Considering the floor operator $\lfloor . \rfloor$ introduced earlier and the total amount of chartists N, the amount of individual agents initially invested in the coupon and dividend paying assets is hence given by:

$$\begin{aligned} N^i_0 = \left\{ \begin{array}{lll} 1 - \sum _{k=2}^{M}{N^{k}_0} &{} \text {for}\quad i = 1, \\ \left\lfloor N \frac{1 - x^{c,M}_0}{M-1} \right\rfloor &{} \text {for}\quad i \in \{2,\ldots ,M-1\}, \\ \Big \lfloor N x^{c,i}_0 \Big \rfloor &{} \text {for}\quad i = M.\end{array}\right. \end{aligned}$$

(76)

This further enables to initialize the opinion indices involved in the allocation process governing the aggregated chartist’s investments:

$$\begin{aligned} s^{\text {step1}}_0= & {} 1 - 2 x^{c,M}_0, \end{aligned}$$

(77)

$$\begin{aligned} s^{\text {step2}}_0 |_{i,i+1}= & {} \frac{1}{2} \frac{N^{i+1}_0 - N^i_0}{N^{i+1}_0 + N^i_0}, \end{aligned}$$

(78)

$$\begin{aligned} s^{\text {step3}}_0 |_{i}= & {} \frac{(M-1) N^i_0 - \sum _{k=1, k \ne i}^{M-1}{N^k_0}}{(M-1)N^i_0 + \sum _{k=1, k \ne i}^{M-1}{N^k_0}}. \end{aligned}$$

(79)

Table 3 Set of constant parameters and initialization of the variables involved in the market model

Table 3 lists the parameters and initial values of the variables involved in the market model. The values included in parentheses, respectively, correspond to each of the assets in terms of increasing order of maturities. The initialization of the herding propensity is achieved as in Kaizoji et al. (2015).

The annualized yields to maturity are successively initialized as follows:

$$\begin{aligned} \text {YTM}^i_0 = \left\{ \begin{array}{lll} \frac{2 \times 250 \times r_f}{2 + 250 \times r_f} &{} \text {for}\quad i = 1 \\ \frac{d^i_0 \times 250 + \frac{\text {FV}^i - P^i_0}{N}}{\frac{P^i_0 + \text {FV}^i}{2}} &{} \text {for}\quad i \in \{2,\ldots ,M-1\}, \\ \frac{d^M_0 \times 250 }{P^M_0} &{} \text {for}\quad i = M.\end{array}\right. \end{aligned}$$

(80)

One can hence observe that the initial yield curve is flat as each yield initially equals $2.5\%$, with the exception of the yield of the perpetual bond initialized as $\text {YTM}^M_0 = 3 \%$. This precise set-up is justified by the need to analyse the impacts of the agents’ investments on a neutral configuration of the yield curve. This further explains the selection of the face values of the corresponding assets.

Building on these developments, the evaluations of the initial artificial changes of the market values, modified and Macaulay durations, convexities, discount factors and implied spot rates of each coupon-paying bond are straightforward. The according values are presented in Table 4 of Appendix A.4.

Finally, given the previous developments, the expectations of the future returns are initialized as:

$$\begin{aligned} \mu ^i_0 = \left\{ \begin{array}{ll} \frac{1}{2} \Big \{ \Big [(1 +\text {YTM}^i_0)^{\frac{1}{250}}-1 \Big ] + \Big [(1 + \text {ISR}^i_0)^{\frac{1}{250}}-1 \Big ] \Big \} + \frac{d^i_0}{P^i_0} &{} \text {for}\quad i \in \{2,\ldots ,M-1\}, \\ \Big [(1 +\text {YTM}^i_0)^{\frac{1}{250}}-1 \Big ] + \frac{d^M_0}{P^M_0} &{} \text {for}\quad i = M.\end{array}\right. \end{aligned}$$

(81)

This concludes the presentation of the initialization of the constant parameters and variables defined previously and leads to the subsequent analyses of the phenomena emerging from numerical simulations of the market model.

3.2 Time series analyses

The following section specifies the dynamics emerging in typical simulations. Figures 2 and 3 show representative time series of the endogenous prices $P^i_t$, yields to maturity $\text {YTM}^i_t$, wealth fractions $x^{f,i}_t$ and $x^{c,i}_t$, transition probabilities $p^{\pm }_t(s^{\text {step}1}_t,H^{\text {step1}}_t)$, momentum $H^{\text {step1}}_t$, and three representative opinion indices $s^{\text {step2}}_t |_{i,i+1}$ created with a constant and a mean-reverting herding propensity. The herding propensity $\kappa _t$ corresponds to the inverse temperature of the underlying Ising-structure of the noise traders’ decision process. In the simulations involving a constant $\kappa$ and where $\kappa < \kappa _c$, the noise traders are in the disordered regime. The stochastic mean reverting herding propensity has the same mean value $\kappa _{\mu }=0.98 \times \kappa _c$ as the constant one, but fluctuates transiently above the critical value. This results in a polarization of the noise traders’ opinions and can be observed through bubbles in the corresponding prices’ time series. These two different types of simulations are analysed hereafter.

The time series of the asset prices $P^i_t$ show the divergence of the preferred stock’s prices from those of the bonds. This divergence is governed by the growth rate of $r_d$. The subcritical herding behaviour associated with the constant nature of the herding propensity results solely in subtle deviations from the average trajectory given by the previous rate. In contrast, the simulations for which the herding propensity varies are found to generate transient regimes when the preferred stock prices undergo turbulent patterns. Such periods are, e.g. observed for $t \in [1300,1700]$. Analyses of the super-exponential growth characterizing these bubbles are provided in Kaizoji et al. (2015), among others. Besides their prices, the bonds are characterized by their yields to maturity $\text {YTM}^i_t$. In both simulations, the yields fluctuate around their initial values with an amplitude related to the maturities of the bonds. These fluctuations reach higher values for bonds of low maturities. This constitutes an emerging phenomenon referred to as the humped “term structure of volatility”. Moreover, according to the well-known price–yield relationship, the yields of the perpetual bond reach their lowest values during the peak of the price bubbles, e.g. for $t \in [1400,1700]$ in Fig. 3. The reverting of the prices of $A^M$ to levels directed by the growth rate of $r_d$ as shown $t \in [3000,3100]$ symmetrically drives the yields back to their initial value of $3 \%$.

Fig. 2

Time dependence of the variables associated with the market model. The simulation length is set to $T_{sim} = 5000$ time steps and the herding propensity remains constant throughout the simulation

Fig. 3

Time dependence of the variables associated with the market model. The simulation length is set to $T_{sim} = 5000$ time steps and the herding propensity follows a stochastic mean-reverting process

The noise traders switch their investments between the different assets with a probability expressed in terms of the momentum of the prices of the preferred stock and of the opinion of others. The dynamics of the momentum $H^M_t$ of the preferred stock, of the probabilities $p^{\pm }_t(s^{\text {step}1}_t,H^{\text {step1}}_t)$ governing the rebalancing between the stock and the fixed income portfolio as well as of three opinion indices $s^{\text {step2}}_t |_{i,i+1}$ are detailed in the following. Overall, the momentum parameter $H^{\text {step1}}_t$ evaluated in terms of the two sub-parameters $H^M_t$ and $H^{\text {step1}}_t$ replicates the main trends of the preferred stock’s prices. This is due to the high amplitudes of the exponentially weighted moving average of the returns of $A^M$ in comparison with the stable spreads proposed by the yield curve during the simulations. The latter amplitudes are nevertheless found to be of approximately an order of magnitude lower than those associated with the opinion index $s^{\text {step1}}_t$. This explains the strong in influence given to the latter variable in the outcome of the Bernoulli trials. The transition probabilities fluctuate symmetrically around $(p^+ + p^-)/2$ and are bounded in $[0,\approx 0.15]$ for the simulation governed by a constant herding propensity. The smooth fluctuations observed in this simulation further explain the lack of abrupt changes in the preferred stock’s price dynamics, as can be expected from the outcome of the Bernoulli trials involving these probabilities and generating no significant rebalancing of the aggregated chartist’s portfolio between the two pools of assets. In contrast, the transition probabilities observed in the simulation involving the stochastic herding propensity increase sharply during the emergence or burst of the bubbles. The noise traders’ opinion formation is illustrated by the following pairs of assets selected to avoid any redundancy in the analyses and to assess the integration of $A^1$ in ${\mathcal {P}}$: ($A^1, A^2$), ($A^7$,$A^1$) and ($A^3$,$A^4$). All three indices are found to fluctuate around a long-run average of zero with the same amplitudes, indicating an equilibrated dispersion of the agents in ${\mathcal {P}}$. The risk-free asset is hence not introducing any bias in the allocation process. This statement is emphasized by the fact that the trends in the fluctuations of the indices $s^{\text {step2}}_t |_{1,2}$ and $s^{\text {step2}}_t |_{7,1}$ replicate each other and are approximately symmetrically mirrored by those of $s^{\text {step2}}_t |_{3,4}$. The emphasis is furthermore directed towards the ascertainment of the abrupt saturation of the opinion indices observed during periods of herding towards the preferred stock proposed by the simulation involving a varying herding propensity. Considering the fact that a bad performing bond is penalized twice: the first time by the loss of agents transiting to the neighbouring bonds and the second time by the outflow of agents towards $A^M$, it appears that a very small amount of agents remain invested in the bonds during the culmination of the preferred stock’s price bubbles. As a result, the opinion indices become very sensitive to the transition of agents as dictated by the outcome of the Bernoulli trials, explaining their abrupt variations at the corresponding time steps.

The transition probabilities define the aggregated noise trader’s investment decisions as they determinate the update of the number of traders invested each asset. The resulting time series of the wealth fraction $x^c_t$ invested in each asset are analysed subsequently. A particular mention tackles the transient dynamics of $x^{c,M}_t$ occurring in the simulation involving the varying $\kappa _t$. For $t \in [2200,2400]$, $s^{\text {step1}}_t$ is found to be entirely polarized at the culmination of the corresponding price bubble. The lock-in effect occurring at this moment is straightforwardly accompanied by the following observations: $x^{c,M}_t = 1$ and $p^{+}_t(s^{\text {step}1}_t,H^{\text {step1}}_t) = 0$. In this case, the reaching of a threshold value by the herding propensity satisfies the statement proposed in Sornette (1994) and granting the responsibility of the “sweeping of an instability” of the latter parameter to instigate a self-reinforcing loop leading to a gradual feedback between $H^M_t$ and $s^{\text {step1}}_t$ provoking an outflow of agents from ${\mathcal {P}}$ to $A^M$. At the apogee of this phenomenon, no agent is spurred to transit to the fixed income portfolio. This situation is, however, unsustainable and the flux of agents between the two asset classes is reversed as soon as $p^+_t > 0$, leading to the inevitable burst of the bubble. The market model is moreover found to be subjected to an initial phase of auto-regulation occurring at the beginning of each simulation. The quintessence of this observation lies in the increase of $x^{c,M}$ from $30 \%$ to a stable level of $\approx 50 \%$ on the long-run. As a consequence of the resulting increase of the preferred stock’s price level, $x^{f,M}$ is further reduced at the benefit of $x^{f,1}$. This is due to the counter-cyclical reaction of the fundamentalist trader to the evolution of the market. More generally, this “price taker” characteristic is found to be revealed in the time series of the corresponding wealth fractions mirroring the fluctuations of the prices of $A^M$.

Figures 4 and 5 provide a comparison of the fundamentalist’s and aggregated chartist’s portfolio statistics for simulations, respectively, involving a constant and a mean-reverting herding propensity. They show the evolution of the average portfolio yields $\text {AvgYieldMV}^i_t$ and $\text {AvgYieldBPV}^i_t$, aggregate wealth levels $W_t^i$ as well as approximated portfolio modified durations $\text {ApproxModDur}^i_t$, convexities $\text {ApproxConv}^i_t$ and dispersions $\text {AvgDisp}^i_t$.

Fig. 4

Time dependence of the variables associated with the market model. The simulation length is set to $T_{sim} = 5000$ time steps and the herding propensity remains constant throughout the simulation

Fig. 5

Time dependence of the variables associated with the market model. The simulation length is set to $T_{sim} = 5000$ time steps and the herding propensity follows a stochastic mean-reverting process throughout the simulation

First, the average market-value-weighted and basis-point-value-weighted portfolio yields to maturity ($\text {AvgYieldMV}^i_t$ and $\text {AvgYieldBPV}^i_t$) are analysed. The portfolio yields computed from the two different approaches provide different assessments of the attractivity of each investment strategy. The MV-weighted yields tend to replicate the patterns of the yields of the preferred stock less than is achieved by the BPV-weighted yields. This is striking for the fundamentalist trader who constrains her endowments to $A^M$ during the simulations. The BPV-weighted yields contrast with the other yields in accordance to the expression of Eq. (113). In particular, as the perpetual bond is granted a much higher duration than the other bonds, its prices are much more sensitive to variations of its yields, explaining why $A^M$ takes a significantly greater share in the sum of the absolute price variations in the determination of the yields of concern.

The average portfolio modified durations $\text {AvgModDur}^i_t$ tend to replicate the dynamics of the modified durations associated with the traders’ predominant investment. The fundamentalist is hence found to have an almost stable portfolio modified duration over time with only few exceptions arising during the transient market unfolding. As a result, the latter variable is hence fluctuating between 20 and 30 years during periods of quietude. The portfolio modified duration of the fundamentalist is also more sensitive to increases of the perpetual bond’s duration than to the according decrease of $x^{f,M}_t$. This is illustrated by the humps in the time series of the corresponding variable generated in the simulation involving the varying herding propensity. The present market model moreover provides a striking observation during periods of extreme transient phenomena and unveiling the limitation of the well-known fixed income pricing formula expressed in Eq. (102). In fact, the relationship implying, among others, that simultaneous increases of the duration and decreases of the yields lead to an increase in the corresponding asset’s market value is found not to hold in turbulent periods for the portfolio of the fundamentalist trader. As one can observe during any of the price bubbles associated with $A^M$, the wealth of the fundamentalist decreases despite the fact that the portfolio modified duration increases and that the portfolio yield decreases. The convexity is besides found to increase, adding no interference to the previous statement. This sheds lights into the boundedness of the earlier mentioned pricing formula and explains why one might loose commonly used landmarks when managing a fundamentalist portfolio during periods of market unrest. Consequently, surfing partly on the price bubble might help the fundamentalist to contain the losses caused by the massive exodus of the chartists from the other assets to the perpetual bond subjected to the bubbles of concern. Following the tendencies observed for the portfolio modified durations, the average portfolio convexities $\text {AvgConv}^i_t$ and dispersions $\text {AvgDisp}^i_t$ propose patterns strikingly sensitive to the dynamics of the perpetual bond’s prices. One can again observe the fact that the rebalancing executed by the fundamentalist as a reaction to the price deviations of $A^M$ is not sufficient to counteract the overall tendencies of the parameters of concern to follow the dynamics of the latter asset.

The comparison of the normalized wealth levels $W^f_t$ of the fundamentalist and $W^c_t$ of the aggregated chartist reveals the fact that the fundamentalist’s strategy outperforms the chartist’s one irrespectively of the process governing the evolution of the herding propensity. The notion of opportunity costs is however crucial to evaluate the reversal of the latter trends on the short-term during turbulent periods. As a matter of fact, such costs appear to be incurred by the fundamentalist during the build-up of the bubbles of the prices of $A^M$. Nevertheless, the subsequent bursts instigated by the flee of individual chartists towards ${\mathcal {P}}$ bringing the variable $x^{c,M}_t$ to anterior levels happen to cancel all the gains of the aggregated chartist, while the fundamentalist kept benefiting from the constant price return of the risk-free asset. As a result, the decreases of $W^f$ observed during the build-up of the price bubbles are adjusted by the still invariable price returns granted by $A^1$ left unimpaired by the variations of $x^{f,1}_t$ and $x^{c,1}_t$. This explains the slight increase of the difference of the wealth levels of the aggregated traders after each burst of a bubble.

4 Analyses

4.1 Application of the Vasicek model of interest rates

The previous developments introduced a constant risk-free rate provided by the zero-coupon bond of one-year maturity. As a matter of fact, the rates provided by such risk-free assets are not constant in reality. The solution of the Vasicek model of interest rates developed by the eponym author in Vasicek (1977) is applied in the market model to provide a varying exogenous arbitrage-free risk-free rate. As a result, the comparison of the emerging dynamics associated with this new set-up with those of the one considering $r_f = \text {cst.}$ shall demonstrate the capacity of the new model to reproduce strikingly similar autocorrelations of the volatilities of the yields to maturity in comparison with those observed in the US Treasury market. The statistical estimation of the parameters involved in the expression of the varying risk-free rate is first achieved below, before the presentation of the emerging dynamics arising in typical simulations. The analyses of the autocorrelations are presented in Sect. 4.2.2.

The “one time step” rates are determined from the following closed-form expression:

$$\begin{aligned} r_{f,t} = r_{f,t-1} e^{- \lambda } + \theta (1 - e^{-\lambda }) + \nu _t \sqrt{\frac{\sigma _r^2 (1 - e^{-2\lambda })}{2 \lambda }}, \end{aligned}$$

(82)

where $\nu _t \sim {\mathcal {N}}(0,1)$ is an i.i.d random variable and $\lambda$, $\theta$ and $\sigma _r$ are estimated in the same way as achieved for the herding propensity. The long-run average of the risk-free rate is equal to the constant value of the rate applied earlier, namely with $\theta = 8 \times 10^{-5}$. The mean-reversion strength is subsequently estimated according to:

$$\begin{aligned} \lambda = \frac{1}{\Delta T} \log _{10} \Big ( \frac{r_{\mathrm{max}} - \theta }{r_{critic} - \theta } \Big ), \end{aligned}$$

(83)

where $r_{\mathrm{max}} = 1 \times 10^{-4}$ and $r_{\mathrm{critic}} = 9.5 \times 10^{-5}$. The value attributed to $\theta$ is known and $\Delta T = 60$ time steps. One can hence evaluate the mean-reversion strength as $\lambda \approx 2 \times 10^{-3}$. Finally, $\sigma _r$ is obtained as a function of the mean-reversion rate as

$$\begin{aligned} \sigma _r = 2 \lambda {\mathrm {Var}}_{r_{f,t}}, \end{aligned}$$

(84)

where the variance of the risk-free rate is such that ${\mathrm {Var}}_{{\mathrm{r}}_{\mathrm{t}}} = 2 \times 10^{-5}$, twice as low as the variance of the perpetual bond’s returns. The initial value of the risk-free rate is finally set to $r_{f,0} = \theta$.

Figure 6 presents the time series generated by a typical simulation involving a varying herding propensity as well as the stochastic process governing the evolution of $r_{f,t}$. The evolution of the risk-free rate has a strong reflection on $x^{f,1}$, while the preferred stock’s prices show a clear symmetric pattern of its dynamics. The evolution of the wealth levels further illustrates the indirect impact of the varying risk-free rate. The reason for the high sensitivity of the fundamentalist’s wealth fractions to the variations of $r_{f,t}$ is straightforward when considering Eqs. (30) and (31). In this regard, the fundamentalists’ investments are hence found to have a major impact on the assets’ prices, being mainly driven by the evolution of the risk-free rate and having a significant impact on the wealth levels of both aggregated traders. There is, however, no form of self-reinforcing feedback arising between the fundamentalist’s investments and the corresponding impacts on the assets’ prices, the latter agent having no incentive to follow series of positive returns and being predominantly reacting to the stochastic variations of the risk-free rate. This further explains why one cannot find clear trends in the evolution of the difference of both agent’s normalized wealth levels, especially after the periods of extreme transient phenomena, where the earlier mention about the opportunity costs is not holding anymore as one is not guaranteed to obtain a stable return in the money market instrument during the culmination of the preferred stock’s price bubbles.

Fig. 6

Time dependence of the variables associated with the market model. $T_{sim} = 5000$ time steps and the herding propensity follows a stochastic mean-reverting process. The evolution of the risk-free rate $r_{f,t}$ is governed by the Vasicek model of interest rates

4.2 Dynamics of the fixed income parameters

The following developments focus on the emerging dynamics of the yields to maturity, implied spot rates and discount factors. For the sake of clarity, the l.h.s. (respectively, r.h.s.) of each illustration is associated with the simulation involving the constant herding propensity (respectively, varying).

Fig. 7

Surface generated from the time dependencies of the yields to maturity of the bonds of finite maturity. The yields of the perpetual bond are identified separately. The l.h.s. (respectively, r.h.s.) is associated with the simulation involving a constant herding propensity (respectively, mean-reverting)

4.2.1 Time dependencies of the fixed income parameters

Figures 7, 8 and 9 present the surfaces generated from the time series associated with the yields to maturity, implied spot rates and discount factors obtained from the two typical simulations presented in Sect. 3.2. The yields of the perpetual bond are represented separately in Fig. 7. Their distinction from the other yields is achieved for the sake of clarity as the corresponding asset is omitted in the computation of the yield curves. An ensemble of $T_{sim} = 5000$ yield curves is generated in each simulation, i.e. one yield curve per time step. The impacts of the aggregated agents’ investments are straightforward considering the initial flat and neutral configuration of the yield curve. The increasing volatile character of the yields to maturity in terms of decreasing maturities is striking in both simulations and further reveals the aptitude of the present market model to produce a humped term structure of volatility. This paradigm is subsequently tackled in Sect. 4.2.2.

Fig. 8

Surface of the implied spot rates of the bonds of finite maturity. The l.h.s. (respectively, r.h.s.) is obtained from a constant herding propensity (respectively, mean-reverting)

Fig. 9

Surface of the discount factors of the bonds of finite maturity. The l.h.s. (respectively, r.h.s.) is obtained from a constant herding propensity (respectively, mean-reverting)

The implied spot rates defined according to Eq. (97) provide accurate approximations of the yields to maturity. As a result, the surfaces presented in Fig. 8 are very similar to those of Fig. 7. The spikes appearing at $t \in [2300,2600]$ in the former very accurately replicate those emerging at the same period in the latter illustration. Figure 9 finally shows the surfaces generated by the linear interpolations of the discount factors associated with each of the bonds of fixed maturities. As one can observe, the disruptions caused by the aggregated chartist’s investments on the assets’ prices have non-negligible effects on these factors for bonds of high maturities.

This concludes the presentation of the evolution of the yields to maturity, implied spot rates and discount factors.

4.2.2 Volatilities of the yields to maturity

Figure 10 presents the surfaces of volatilities of the yields to maturity. The volatilities are computed over moving windows of 125 time steps in both simulations. They are evaluated at each time step from the following definition,

$$\begin{aligned} \sigma _{\text {YTM}^i,t} = \left. \sqrt{\big \langle (\text {YTM}^i_k - \langle \text {YTM}^i_k \rangle )^2 \big \rangle } \right| ^{t}_{k = t - \tau }, \end{aligned}$$

(85)

where the sample average operator is represented by $\langle . \rangle$ and the selected moving window by $\tau = 125$.

Figure 11 shows the realized volatilities of the yields of US Treasury bonds evaluated between the 18th of November 1993 and the 6th of June 2020, considering moving windows of 250 business days. The selection of this moving window is justified by the need to amplify the trends of the fluctuations of the realized volatilities. Notwithstanding the fact that the volatilities of the yields of the US Treasury bonds are subjected to the exogenous policies of the Federal Reserve, the comparison between those realized endogenously in the present model with those of the US Treasury market enables to assess the capacity of the present market model to account for the main trends observed in reality.

Fig. 10

Surface of volatilities of the yields to maturity in our ABM. The l.h.s. (respectively, r.h.s.) is obtained from a constant $\kappa _t$ (respectively, mean-reverting)

Fig. 11

Surface of volatilities of the yields to maturity of the corresponding US Treasury bonds realized between the 18nth of November 1993 and the 6th of June 2020

The model produces a humped term structure of volatility in both simulations. The term structure of volatility is further decreasing when considering the bonds of maturities larger than one year. Such phenomena are found episodically in the US Treasury market, e.g. between 2004 and 2008, corresponding to $t \in [3200,3950]$ in Fig. 11. On the other hand, the exogenous imposition of low discount rates by the Federal Reserve are commonly known to be reflected on the yields of the US Treasury bonds of short maturities, explaining the boundedness of their volatilities during the corresponding time frames ($t \in [4000,5500]$). The increases of these bounds hence result in increasing volatilities of the yields of bonds of short maturities. As revealed in Cieslak and Povala (2016) and Bertocchi et al. (2005), the volatilities of the US Treasury yield curve are found to be following such a hump-like shape on average, comforting the assertion of the capacity of the present market model to account for the main trends observed for the volatilities of the yields of US Treasury bonds.

4.2.3 Autocorrelations of the volatilities of the yields to maturity

Figures 12 and 13 present the surfaces generated from the coefficients of autocorrelation of the volatilities of the yields to maturity computed for each bond, excluding the risk-free one and including the perpetual one. The former and the latter, respectively, consider the simulations involving a constant and a varying risk-free rate. The surfaces are generated by computing the average of the coefficients of autocorrelations from a set of 300 simulations. The first 500 steps of each simulation are furthermore omitted in the computations due to the auto-regulation phenomenon commented in Sect. 3.2. The above-mentioned coefficients are determined for each simulation from Eq. (86), considering the volatilities $\sigma _{\text {YTM}^i,t}$ determined earlier and the lags $l \in [0,500]$:

$$\begin{aligned} \begin{aligned} \text {ACF}_l\left( \sigma _{\text {YTM}^i}\right)&= \frac{{{\mathrm {Cov}}}\left[ \sigma _{\text {YTM}^i,t}\sigma _{\text {YTM}^i,t-l}\right] }{\sqrt{{\mathrm {Var}}\left[ \sigma _{\text {YTM}^{\mathrm{i}},{\mathrm{t}}}\right] {\mathrm {Var}}\left[ \sigma _{\text {YTM}^{\mathrm{i}},{\mathrm{t-l}}}\right] }}, \\&=\frac{\left\langle \sigma _{\text {YTM}^i,t}\sigma _{\text {YTM}^i,t-l} \right\rangle - \left\langle {\sigma ^2_{\text {YTM}^i,t}} \right\rangle }{\left\langle \sigma ^2_{\text {YTM}^i,t} \right\rangle - {\left\langle \sigma _{\text {YTM}^i,t} \right\rangle }^2 }. \end{aligned} \end{aligned}$$

(86)

Fig. 12

Surface of autocorrelation associated with the simulation involving a constant risk-free rate and a constant $\kappa _t$ on the l.h.s. (respectively, mean-reverting $\kappa _t$ on the r.h.s.)

Fig. 13

Surface of autocorrelation associated with the simulation involving a mean-reverting risk-free rate and a constant $\kappa _t$ on the l.h.s. (respectively, mean-reverting $\kappa _t$ on the r.h.s.)

The comparison between Figs. 12 and 13 reveals the fact that the surfaces associated with the market model involving a varying risk-free rate decay much faster than the others in terms of increasing lags. This emerging phenomenon is hence attributed to the stochastic character of the evolution of the risk-free rate. There is moreover no difference on the aggregated level for different selections of the processes governing the evolution of the herding propensity.

The comparison of the surfaces presented in Figs. 12 and 13 with the surface of autocorrelation associated with the US Treasury market provided in Fig. 14 and computed from the earlier-mentioned volatilities further enables to distinguish the capacity of the simulations involving the varying risk-free rates to replicate the same emerging phenomena as they reproduce similar autocorrelation coefficients at the corresponding lags. More precisely, the latter simulations are capable to replicate the steeper gradients of the autocorrelation surfaces until the 250th lag as well as the milder ones for higher lags directing the surfaces towards null autocorrelations at the 500th lag.

Fig. 14

Surface of autocorrelation computed on the volatilities of the realized US Treasury yields between the 18nth of November 1993 and the 6th of June 2020

Equation 86 defining the autocorrelations of the volatilities of the yields to maturity aims to uncover persistent dynamics revealing a self-dependency of the observed variable. For example, high autocorrelation values indicate that periods of high volatility induce high volatilities at a later point of time. As one can expect, letting the risk-free rate follow an independent stochastic process affects the dynamics of the yields to maturity through the transmission channel of the noise trader’s investments presented in step 2 of Sect. 2.4. As a result, the residual self-dependencies of the aforementioned variable are found to diminish at a higher pace when increasing the lags, as one can expect the autocorrelations to be dampened. This is consistent with the empirical observations provided in Fig. 14. However, the latter illustration reveals the fact that the yield volatilities of the bonds of shorter maturities exhibit larger auto-correlations than those of longer maturities. This stylized fact is not reproduced by the present model, as one would need to consider the classical theories of the term structure of interest rates in the traders’ investment processes. According to these theories, fixed income investors choose maturities based on their expectations of future short-term interest rates in order to maximize their overall horizon yield. The relative merits of expectations theory, segmented market theory or liquidity preference theory remains to be decided, as none of these theories have yet gained the upper hand, as exposed in Smith (2014).

5 Conclusion

The agent-based market model proposed in this research has extended the model proposed in Kaizoji et al. (2015) by introducing several bonds of fixed maturities. The inclusion of the various fixed income parameters allowed us to describe the investment processes of both aggregated traders to the peculiarities of the market in which they evolve. The market set-up initialized by a flat yield curve demonstrated its capacity to reproduce several stylized facts of the US Treasury market. The simulations involving a varying herding propensity showed that the well-known bond-pricing formula involving the durations, convexities and shifts in yields to maturity does not hold at the fundamentalist’s portfolio level during periods of market unrest. Subsequently, the implementation of Vasicek’s model of interest rates to simulate the evolution of the risk-free rate was found to allow for the linearly dependent fundamentalist’s wealth fractions to replicate its dynamics and have a significant impact on the asset prices. The analyses of the volatilities of the yields to maturity further enabled to assess the capacity of the market model to produce a humped term structure of volatility as is observed on average in the US Treasury market. Finally, the market model introducing the mean-reverting process governing the dynamics of the risk-free rate also proved superior to the one involving $r_f = \text {cst.}$ with regards to its capacity to reproduce the main characteristics of the surface of autocorrelation of the volatilities of the yields to maturity of the US Treasury bonds for selected the time-frame.

Appendix

1.1 Parameters of the fixed income market model

• Security level

The following developments introduce the parameters mentioned in Sect. 2.1. First and foremost, the famous relationship between an asset’s no-arbitrage price ${\text {PV}}$, face value ${\text {FV}}$, coupon payments $\text {C}$, time to maturity N and yield to maturity ${\text {YTM}}$ is expressed as:

(87)

Representing the discount factor associated with the payments occurring at maturity by $\delta (N)$, the no-arbitrage price can be further expressed as:

(88)

exacerbating the expression of the bonds’ dirty prices in terms of those of a perpetuity and of a zero-coupon bond as follows:

$$\begin{aligned} {\text {PV}} = {\text {PV}}_{\text {perpetuity}} (1 - \delta (N)) + \delta (N) {\text {PV}}_{\text {zero-coupon bond}}. \end{aligned}$$

(89)

Yields to maturity Considering the fact that each asset $A^i$ for $i \in \{2,\ldots ,M\}$ pays a coupon or a dividend of $d_t$ at each time step and consisting of the accrued interest associated with $\mathrm{d}t = 1$, the corresponding yields to maturity are evaluated as:

$$\begin{aligned} \text {YTM}^i_t \approx \frac{250 \times d^i_t + \left( \frac{\text {FV}^i - {\text {PV}}^i_t}{N^i}\right) }{\frac{\text {FV}^i + {\text {PV}}^i_t}{2}}, \end{aligned}$$

(90)

where the day-count convention defines a year as consisting of 250 trading days. The accrued interests $\text {AI}^i_t$ of each asset $A^i$ are further determined at each time step t by:

$$\begin{aligned} \text {AI}^i_t = \text {FV}^i \times \frac{\text {C}^i}{\text {FV}^i} \times \frac{\mathrm{d}t}{T_{\text {year}}}, \end{aligned}$$

(91)

where $T_{\text {year}} = 250$. Accordingly, the running yield of the preferred stock $A^M$ is given by:

$$\begin{aligned} \text {YTM}^M_t \approx \frac{2 \times 250 \times d^M_t}{\text {FV}^M+{\text {PV}}^M_t}. \end{aligned}$$

(92)

Discount factors The discount factors provide a practical bootstrapping technique upon which the spot rates are determined. The methodology applied to determine the discount factors is detailed hereafter:

• 1. Start by determining the discount factor of the risk-free bond:

  $$\begin{aligned} \text {DF}^1_t = \frac{{\text {PV}}^1_t}{\text {FV}^1}. \end{aligned}$$

  (93)

  Setting ${\text {PV}}^1 = 1$ and $\text {FV}^1 = 1 + 250 \times r_f$, one obtains:

  $$\begin{aligned} \text {DF}^1_t = \frac{1}{1 + 250 \times r_f} \quad = \text {cst.} \quad \forall t. \end{aligned}$$

  (94)

• 2. Determinate the remaining $\text {DF}^i_t$ in an increasing order of maturities:

$$\begin{aligned} \text {DF}^i_t = \frac{{\text {PV}}^i_{t} - 250 \times d^i_{t} \sum _{k=1}^{i-1}{\text {DF}^k_t}}{250 \times d^i_{t} + \text {FV}^i}. \end{aligned}$$

(95)

As one can observe from Eq. (95), discount factors can only be assigned to bonds of fixed maturities. These factors are further related to the implied spot rates $\text {ISR}^i_t$ by:

$$\begin{aligned} \text {DF}^i_t = \frac{1}{\Big (1 + \frac{\text {ISR}^i_t}{\text {PER}} \Big )^{N \times \text {PER}}}, \end{aligned}$$

(96)

where the periodicity accounting for semiannual coupon payments is set to $\text {PER} = 2$.

Implied spot rates The implied spot rates are accordingly isolated by the subsequent expression:

$$\begin{aligned} \text {ISR}^i_t = \left[ \left( \frac{1}{\text {DF}^i_t} \right) ^{\frac{1}{N^i \times \text {PER}}} -1 \right] \times \text {PER}. \end{aligned}$$

(97)

Variations of market values Assuming that the prices of each bond are continuous and twice differentiable w.r.t. to the yields to maturity, we have

$$\begin{aligned} \Delta \text {MV}^i_t \approx \frac{\partial \text {MV}^i_t}{\partial \text {YTM}^i_t} d \text {YTM}^i_t + \frac{1}{2} \frac{\partial ^2 \text {MV}^i_t }{\partial {\text {YTM}^i_t}^2} (d \text {YTM}^i_t)^2, \end{aligned}$$

(98)

from which the following parameters are introduced, considering only flat variations of the yield curve.

Durations and convexity

The Macaulay duration is first expressed as:

$$\begin{aligned} \text {MacDur}^i_t = - \frac{\partial \text {MV}^i_t}{\partial \text {YTM}^i_t} \frac{1 + \text {YTM}^i_t}{\text {MV}^i_t}, \end{aligned}$$

(99)

followed by the definition of the modified duration:

$$\begin{aligned} \text {ModDur}^i_t = - \frac{\partial \text {MV}^i_t}{\partial \text {YTM}^i_t} \frac{1 }{\text {MV}^i_t} = \frac{\text {MacDur}^i_t}{1 + \text {YTM}^i_t}. \end{aligned}$$

(100)

The convexity is finally formulated as:

$$\begin{aligned} \text {Conv}^i_t = \frac{\partial ^2 \text {MV}^i_t}{\partial {\text {YTM}^i_t}^2} \frac{1}{\text {MV}^i_t}, \end{aligned}$$

(101)

allowing to express Eq. (98) in terms of the above-introduced parameters:

$$\begin{aligned} \frac{d \text {MV}^i_t}{\text {MV}^i_t} \approx -(\text {ModDur}^i_t \times \mathrm{d}\text {YTM}^i_t ) + \frac{1}{2}(\text {Conv}^i_t \times (\mathrm{d}\text {YTM}^i_t)^2). \end{aligned}$$

(102)

These parameters are further evaluated from the approximations provided in Smith (2014) and are first applied to individual securities before being associated with the aggregated traders’ portfolios.

The modified duration is approximated by:

$$\begin{aligned} \text {ApproxModDur}^i_t = \frac{\Delta \text {MV}^i_t |_{\text {YTM}^i_t \downarrow } - \Delta \text {MV}^i_t |_{\text {YTM}^i_t \uparrow }}{2 \times \Delta \text {YTM}^i_t \times \text {MV}^i_t}, \end{aligned}$$

(103)

where $\Delta \text {MV}^i_t |_{\text {YTM}^i_t \downarrow (\uparrow )}$ represents the variation of the bond’s market value resulting from an artificial decrease (increase) of the yield by $- (+) 0.001$.

The convexity is furthermore approximated by:

$$\begin{aligned} \text {ApproxConv}^i_t = \frac{\Delta \text {MV}^i_t |_{\text {YTM}^i_t \downarrow } + \Delta \text {MV}^i_t |_{\text {YTM}^i_t \uparrow } - 2 \times \text {MV}^i_t}{ {\Delta \text {YTM}^i_t}^2 \times \text {MV}^i_t}. \end{aligned}$$

(104)

In particular, considering the absence of arbitrage, the assets’ market values and their corresponding variations arising from artificial shifts of the yields are determined according to:

$$\begin{aligned} \text {MV}^i_t= & {} \frac{ 250 \times d^i_t}{\text {YTM}^i_t} \left( 1 - \frac{1}{(1+\text {YTM}^i_t)^N} \right) + \frac{\text {FV}^i}{(1+\text {YTM}^i_t)^N}, \end{aligned}$$

(105)

$$\begin{aligned} \text {MV}^i_t |_{\text {YTM}^i_t \downarrow }= & {} \frac{ 250 \times d^i_t}{\text {YTM}^i_t - 0.001} \left( 1 - \frac{1}{(1+\text {YTM}^i_t - 0.001)^N} \right) \nonumber \\&+ \frac{\text {FV}^i}{(1+\text {YTM}^i_t - 0.001)^N}, \end{aligned}$$

(106)

$$\begin{aligned} \text {MV}^i_t |_{\text {YTM}^i_t \uparrow }= & {} \frac{ 250 \times d^i_t}{\text {YTM}^i_t + 0.001} \left( 1 - \frac{1}{(1+\text {YTM}^i_t + 0.001)^N} \right) \nonumber \\&+ \frac{\text {FV}^i}{(1+\text {YTM}^i_t + 0.001)^N}. \end{aligned}$$

(107)

The Macaulay duration is finally approximated by:

$$\begin{aligned} \text {ApproxMacDur}^i_t = \text {ApproxModDur}^i_t (1 + \text {YTM}^i_t). \end{aligned}$$

(108)

• Portfolio level

Durations and convexity The market-value-weighted (MV-weighted) average modified portfolio durations are successively expressed as:

$$\begin{aligned} \begin{aligned} \text {AvgModDur}^i_t&= \sum _{k=1}^{M}{x^{i,k}_{t} \times \text {ModDur}^{k}_t} \\&\approx \sum _{k=1}^{M}{x^{i,k}_{t} \times \text {ApproxModDur}^k_t}, \end{aligned} \end{aligned}$$

(109)

where $i \in \{f,c\}$ depicts the trader of concern. Similarly, the market-value-weighted average Macaulay portfolio durations are introduced as:

$$\begin{aligned} \begin{aligned} \text {AvgMacDur}^i_t&= \sum _{k=1}^{M}{x^{i,k}_{t} \times \text {MacDur}^{k}_t} \\&\approx \sum _{k=1}^{M}{x^{i,k}_{t} \times \text {ApproxMacDur}^k_t}. \end{aligned} \end{aligned}$$

(110)

The average portfolios convexities are further determined as:

$$\begin{aligned} \begin{aligned} \text {AvgConv}^i_t&= \sum _{k=1}^{M}{x^{i,k}_{t} \times \text {Conv}^k_t} \\&\approx \sum _{k=1}^{M}{x^{i,k}_{t} \times \text {ApproxConv}^k_t}. \end{aligned} \end{aligned}$$

(111)

Yields to maturity The market-value-weighted average portfolio yields are evaluated by:

$$\begin{aligned} \text {AvgYieldMV}^i_t = \sum _{k=1}^{M}{x^{i,k}_{t} \times \text {YTM}^k_t}. \end{aligned}$$

(112)

The basis point values (BPVs) of each asset $A^k$ are determined by each trader i as:

$$\begin{aligned} \begin{aligned} \text {BPV}^{i,k}_t&= \text {ModDur}^k_t \times x^{i,k}_t \times W^i_t \times 0.0001 \\&\approx \text {ApproxModDur}^k_t \times x^{i,k}_t \times W^i_t \times 0.0001. \end{aligned} \end{aligned}$$

(113)

Hence, the average basis point value-weighted (BPV-weighted) yields of each trader are evaluated as:

$$\begin{aligned} \text {AvgYieldBPV}^i_t = \frac{\sum _{k=1}^{M}{\text {BPV}^{i,k}_t \times \text {YTM}^k_t}}{{\sum _{k=1}^{M}{\text {BPV}^{i,k}_t}}} . \end{aligned}$$

(114)

Dispersion Ultimately, the average portfolios dispersions are introduced as:

$$\begin{aligned} \text {AvgDisp}^i_t = \text {AvgConv}^i_t (1+ \text {AvgYieldBPV}^i_t)^2 - {\text {AvgMacDur}^i_t}^2 - \text {AvgMacDur}^i_t. \end{aligned}$$

(115)

1.2 Solution of the generalized mean-variance optimization problem

The stochastic optimal control problem assigned to the fundamentalist trader is tackled in the following text. The optimal feedback solution is derived after the formulation of the well-known Hamilton–Jacobi–Bellman (HJB) equation associated with the problem of concern. The solution is furthermore given by the optimal control trajectory $\widehat{{\mathbf {x}}} :[0,T] \mapsto {\mathbb {R}}^{M-1}$ at the origin of the optimal state trajectory $\widehat{W^f}: [0,T] \mapsto {\mathbb {R}}$. The former is expressed as:

$$\begin{aligned} \widehat{{\mathbf {x}}}(t) = \tilde{{\mathbf {x}}}(\widehat{W^f}(t),\nabla _{W^f} {\mathcal {J}}(\widehat{W^f}(t),t),t) \quad \text {for} \quad t \in [0,T]. \end{aligned}$$

(116)

In a nutshell, the control ${\mathbf {x}} = {\mathbf {x}}(W^f,{\mathcal {J}}_{W^f},{\mathcal {J}}_{W^fW^f},t)$ which maximizes the r.h.s. of the HJB Eq. (117) involving the cost-to-go function ${\mathcal {J}}(W^f,t)$ defined in Eq. (24) is first sought for. The resulting expression of the control law is then inserted back into the HJB equation before the resulting PDE expressed in terms of the cost-to-go function is solved. The expression of ${\mathcal {J}}(W^f,t)$ is finally replaced in the expression of the control parameter ${\mathbf {x}}$, unveiling the formulation of the optimal control law.

The HJB equation associated with the present framework is expressed in Eq. (117). For a detailed derivation of the HJB equation through Bellman’s dynamic programming approach and to better grasp its application to the present framework, the reader is invited to refer to Yong and Zhou (1999) or Fleming and Rishel (1975).

$$\begin{aligned} - {\mathcal {J}}_t = \max _{ {\mathbf {x}} \in {\mathbb {R}}^{M-1}} \left\{ {\mathcal {J}}_{W^f} W^f \Big [r_f + {\mathbf {x}}^\top ({\mu } - {\mathbf {e}}r_f) \Big ]+ \frac{1}{2}{\mathcal {J}}_{W^f W^f}{W^f}^2 {\mathbf {x}}^\top {\sigma } {\sigma }^\top {\mathbf {x}} \right\} . \end{aligned}$$

(117)

Recalling the fact that the covariance matrix associated with the SDE of the wealth update is expressed as ${{\varvec{\Sigma }}} = {\sigma } {\sigma }^\top$, Eq. (117) can be reduced to:

$$\begin{aligned} \frac{d}{d{\mathbf {x}}} \left\{ {\mathcal {J}}_{W^f} W^f \Big [r_f + {\mathbf {x}}^\top ({\mu } - {\mathbf {e}}r_f) \Big ]+ \frac{1}{2} {\mathcal {J}}_{W^f W^f}{W^f}^2 {\mathbf {x}}^\top {{\varvec{\Sigma }}} {\mathbf {x}} \right\} = 0, \end{aligned}$$

(118)

which further yields the following expression:

$$\begin{aligned} ({\mu } - {\mathbf {e}}r_f)W^f {\mathcal {J}}_{W^f} + {\mathcal {J}}_{W^f W^f} {W^f}^2 {{\varvec{\Sigma }}} {\mathbf {x}} = 0. \end{aligned}$$

(119)

Considering the latter Eq. (119), one can obtain the formulation of the optimal control law expressed in terms of the first- and second-order derivatives of the cost-to-go function w.r.t. $W^f$:

$$\begin{aligned} \widehat{{\mathbf {x}}} = - \frac{{\mathcal {J}}_{W^f}}{W^f {\mathcal {J}}_{W^f W^f}} {{{\varvec{\Sigma }}}}^{-1} ({\mu } - {\mathbf {e}}r_f). \end{aligned}$$

(120)

Inserting the expression of Eq. (120) into the HJB Eq. (117), one obtains the following PDE expressed in terms of the cost-to-go function:

$$\begin{aligned} {\mathcal {J}}_t + W^f r_f {\mathcal {J}}_{W^f} - \frac{{{\mathcal {J}}_{W^f}}^2}{2 {\mathcal {J}}_{W^f W^f}} ({\mu } - {\mathbf {e}}r_f)^\top {{\varvec{\Sigma }}}^{-1} ({\mu } - {\mathbf {e}}r_f) = 0. \end{aligned}$$

(121)

The following Ansatz introduces the scalar function $h(t) : {\mathbb {R}} \mapsto {\mathbb {R}}$ subjected to the terminal condition given by $h(T) = 1$:

$$\begin{aligned} {\mathcal {J}}(W^f,t)= h(t) \frac{{W^f}^{1 - \gamma }}{1 - \gamma }. \end{aligned}$$

(122)

The relevant derivatives of the value function are consequently obtained as:

$$\begin{aligned} {\left\{ \begin{array}{ll} {\mathcal {J}}_t(W^f,t)= {\dot{h}}(t) \frac{{W^f}^{1 - \gamma }}{1 - \gamma }, \\ {\mathcal {J}}_{W^f}(W^f,t)= h(t){W^f}^{- \gamma },\\ {\mathcal {J}}_{W^f W^f}(W^f,t)= - \gamma h(t) {W^f}^{-(\gamma +1)}, \end{array}\right. } \end{aligned}$$

(123)

allowing to express the solution of the optimization problem ruling the investments of the fundamentalist trader as:

$$\begin{aligned} \widehat{{\mathbf {x}}} = \frac{1}{\gamma } {{{\varvec{\Sigma }}}}^{-1} ({\mu } - {\mathbf {e}}r_f) \end{aligned}$$

(124)

This expression solves the generalized wealth allocation problems faced by any trader subjected to the necessity to maximize an expected CRRA utility function w.r.t. her investments in $M-1$ coupon and dividend paying assets. The rest is allocated in the risk-free asset according to the boundary condition expressed in Eq. (74).

1.3 Market clearing process of the extended market model

The developments provided below shed lights into the market clearing process achieved at each time step by the Walrasian auctioneer. As has been expressed in Sect. 2.5, Eq. (63) is the condition at the origin of the updates of the assets’ prices. The expressions of the excess demands for each asset $A^i$ with $i \in \{2,\ldots ,M\}$, respectively, formulated by the fundamentalist and chartist traders are:

$$\begin{aligned} \Delta D^{f,i}_{t-1 \rightarrow t}= & {} W^f_{t-1} \left( A^i_t + \frac{B^i_t}{P^i_t} \right) \left[ \left( 1 - \sum _{k=2}^{M}{x^{f,k}_{t-1}} \right) (1+r_f) + \sum _{k=2}^{M}{x^{f,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] \nonumber \\&- x^{f,i}_{t-1} W^f_{t-1} \frac{P^i_t}{P^i_{t-1}}, \end{aligned}$$

(125)

$$\begin{aligned} \Delta D^{c,i}_{t-1 \rightarrow t}= & {} W^c_{t-1} x^{c,i}_{t} \left[ \left( 1 - \sum _{k=2}^{M}{x^{c,k}_{t-1}} \right) (1+r_f) + \sum _{k=2}^{M}{x^{c,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] \nonumber \\&- x^{c,i}_{t-1} W^c_{t-1} \frac{P^i_t}{P^i_{t-1}}, \end{aligned}$$

(126)

where both $A^i_t$ and $B^i_t$ are defined according to Eqs. (57) and (58). For the sake of clarity, let $\left( 1 - \sum _{k=2}^{M}{x^{f,k}_{t-1}} \right) (1+r_f)$ be represented by $\Pi _{t-1}$ and $\Phi _{t-1}$ account for $\left( 1 - \sum _{k=2}^{M}{x^{c,k}_{t-1}} \right) (1+r_f)$. The condition of the Walrasian auctioneer is refined hereafter, considering Eqs. (125) and (126):

$$\begin{aligned}&W^f_{t-1} \left( A^i_t + \frac{B^i_t}{P^i_t} \right) \left[ \Pi _{t-1} + \sum _{k=2}^{M}{x^{f,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] - x^{f,i}_{t-1} W^f_{t-1} \frac{P^i_t}{P^i_{t-1}} \nonumber \\&\quad \quad + W^c_{t-1} x^{c,i}_{t} \left[ \Phi _{t-1} + \sum _{k=2}^{M}{x^{c,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] - x^{c,i}_{t-1} W^c_{t-1} \frac{P^i_t}{P^i_{t-1}} = 0.\nonumber \\ \end{aligned}$$

(127)

One can further factorize the terms accounting for the ratios of the prices evaluated at t and $t-1$ by isolating $P_t$ as follows:

$$\begin{aligned}&W^f_{t-1} \left( A^i_t + \frac{B^i_t}{P^i_t} \right) \left[ \Pi _{t-1} + \sum _{k=2}^{M}{x^{f,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] \nonumber \\&\quad \quad + W^c_{t-1} x^{c,i}_{t} \left[ \Phi _{t-1} + \sum _{k=2}^{M}{x^{c,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] \nonumber \\&\quad \quad - P^i_t \left( \frac{x^{c,i}_{t-1} W^c_{t-1} + x^{f,i}_{t-1} W^f_{t-1}}{P^i_{t-1}} \right) = 0. \end{aligned}$$

(128)

Let the parameter $\beta _i$ associated with asset $A^i$ be introduced as:

$$\begin{aligned} \beta _i = \frac{x^{c,i}_{t-1} W^c_{t-1} + x^{f,i}_{t-1} W^f_{t-1}}{P^i_{t-1}}, \end{aligned}$$

(129)

Equation (128) can be consequently expressed as:

$$\begin{aligned}&W^f_{t-1} \left( P^i_t A^i_t + B^i_t \right) \left[ \Pi _{t-1} + \sum _{k=2}^{M}{x^{f,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] \nonumber \\&\quad \quad + P^i_t W^c_{t-1} x^{c,i}_{t} \left[ \Phi _{t-1} + \sum _{k=2}^{M}{x^{c,k}_{t-1}} \left( \frac{P^k_t + d^k_t}{P^k_{t-1}} \right) \right] \nonumber \\&\quad \quad - {P^i_t}^2 \beta _i = 0, \end{aligned}$$

(130)

before being arranged in the following form:

$$\begin{aligned}&P^i_t \left[ W^f_{t-1} A^i_t \Pi _{t-1} + W^c_{t-1} x^{c,i}_{t} \Phi _{t-1} + W^f_{t-1} A^i_t \sum _{k=2}^{M}{x^{f,k}_{t-1}\left( \frac{d^k_t}{P^k_{t-1}} \right) }\right. \nonumber \\&\quad \quad \left. + W^c_{t-1} x^{c,i}_{t} \sum _{k=2}^{M}{x^{c,k}_{t-1}} \left( \frac{d^k_t}{P^k_{t-1}} \right) + \sum _{k=2}^{M}{ \frac{P^k_t}{P^k_{t-1}} \left( x^{f,k}_{t-1} W^f_{t-1} A^i_t + x^{c,k}_{t-1} W^c_{t-1} x^{c,i}_{t} \right) }\right] \nonumber \\&\qquad + W^f_{t-1} B^i_t \Pi _{t-1} + W^f_{t-1} B^i_t \sum _{k=2}^{M}{x^{f,k}_{t-1}} \left( \frac{P^k_t}{P^k_{t-1}}\right) + W^f_{t-1} B^i_t \sum _{k=2}^{M}{x^{f,k}_{t-1}} \left( \frac{d^k_t}{P^k_{t-1}}\right) - {P^i_t}^2 \beta _i = 0.\nonumber \\ \end{aligned}$$

(131)

The parameter $\alpha _{ij}$ associated with the pair of assets i and j is introduced as follows:

$$\begin{aligned} \alpha _{ij} = \frac{x^{f,j}_{t-1} W^f_{t-1} A^i_t + x^{c,j}_{t-1} W^c_{t-1} x^{c,i}_{t}}{P^j_{t-1}} \end{aligned}$$

(132)

and the parameter $\zeta _i$ is defined such that:

$$\begin{aligned} \zeta _i&= W^f_{t-1} A^i_t \Pi _{t-1} + W^c_{t-1} x^{c,i}_{t} \Phi _{t-1} + W^f_{t-1} A^i_t \sum _{k=2}^{M}{x^{f,k}_{t-1}\left( \frac{d^k_t}{P^k_{t-1}} \right) }\nonumber \\&\quad + W^c_{t-1} x^{c,i}_{t} \sum _{k=2}^{M}{x^{c,k}_{t-1}} \left( \frac{d^k_t}{P^k_{t-1}} \right) \nonumber \\&= W^f_{t-1} A^i_t \left( 1 - \sum _{k=2}^{M}{x^{f,k}_{t-1}} \right) (1+r_f) + W^c_{t-1} x^{c,i}_{t} \left( 1 - \sum _{k=2}^{M}{x^{c,k}_{t-1}} \right) (1+r_f) \nonumber \\&\quad + \sum _{k=2}^{M}{\left[ \frac{d^k_t}{P^k_{t-1}} \left( x^{f,k}_{t-1} W^f_{t-1} A^i_t + x^{c,k}_{t-1} W^c_{t-1} x^{c,i}_{t} \right) \right] }. \end{aligned}$$

(133)

Equation (131) can hence be written in the following compact form:

$$\begin{aligned}&P^i_t \left( \zeta _i + \sum _{k=2}^{M}{ P^k_t \alpha _{ik}} \right) + \sum _{k=2}^{M}{\frac{P^k_t}{P^k_{t-1}} W^f_{t-1} B^i_t x^{f,k}_{t-1}} \nonumber \\&\quad \quad + W^f_{t-1} B^i_t \Pi _{t-1} + W^f_{t-1} B^i_t \sum _{k=2}^{M}{x^{f,k}_{t-1}} \left( \frac{d^k_t}{P^k_{t-1}}\right) \nonumber \\&\quad \quad - {P^i_t}^2 \beta _i = 0. \end{aligned}$$

(134)

Successively, let $\chi _{ij}$ be defined as follows:

$$\begin{aligned} \chi _{ij} = \frac{W^f_{t-1} B^i_t x^{f,k}_{t-1}}{P^k_{t-1}} \end{aligned}$$

(135)

and $\lambda _i$ be introduced by the subsequent expression:

$$\begin{aligned} \lambda _{i}= & {} W^f_{t-1} B^i_t \Pi _{t-1} + W^f_{t-1} B^i_t \sum _{k=2}^{M}{x^{f,k}_{t-1}} \left( \frac{d^k_t}{P^k_{t-1}}\right) \nonumber \\= & {} W^f_{t-1} B^i_t \left[ \left( 1 - \sum _{k=2}^{M}{x^{f,k}_{t-1}} \right) (1+r_f) + \sum _{k=2}^{M}{x^{f,k}_{t-1}} \frac{d^k_t}{P^k_{t-1}} \right] . \end{aligned}$$

(136)

Inserting these latter parameters in Eq. (134) further enables to obtain the following expression:

$$\begin{aligned} P^i_t \zeta _i + P^i_t \sum _{k=2}^{M}{ P^k_t \alpha _{ik}} + \sum _{k=2}^{M}{P^k_t \chi _{ik}} - {P^i_t}^2 \beta _i + \lambda _i = 0. \end{aligned}$$

(137)

By isolating the nonzero terms when inserting the Kronecker delta in the summations, one can finally express Equation (137) as:

$$\begin{aligned} {P^i_t}^2 \left( \alpha _{ii} - \beta _i \right) + P^i_t \left( \zeta _i + \chi _{ii} \right) \sum _{k=2, k \ne i}^{M}{ P^k_t P^i_t \alpha _{ik}} + \sum _{k=2, k \ne i}^{M}{P^k_t \chi _{ik}} + \lambda _i = 0. \end{aligned}$$

(138)

This concludes the presentation of the solvable set of $M-1$ nonlinear equations expressed in terms of the $M-1$ asset prices remaining to be evaluated at each time step of the simulations of the model.

1.4 Initialization of the parameters of the fixed income Market

Table 4 below introduces the initial values of the parameters associated with the fixed income market framework. As one could anticipate, the initial values of the modified and Macaulay durations approximate the maturities of the bonds.

Table 4 Initialization of the parameters of the fixed income market. The values correspond, respectively, to each asset in terms of increasing maturities, starting with the one-year zero-coupon risk-free bond and ending with the perpetual bond