Linear Stochastic Dividend Model

TL;DR

This paper introduces a new polynomial process-based model for stock and dividend derivatives, enabling explicit computation of moments and closed-form pricing of futures, with accurate option pricing via moment matching.

Contribution

It proposes a novel stochastic dividend model that is tractable and non-affine, allowing explicit moment calculations and improved derivative pricing methods.

Findings

Closed-form expressions for stock and dividend futures prices.

Accurate approximation of stock and dividend options using moment matching.

Model maintains positivity of stock and non-negativity of dividends.

Abstract

In this paper we propose a new model for pricing stock and dividend derivatives. We jointly specify dynamics for the stock price and the dividend rate such that the stock price is positive and the dividend rate non-negative. In its simplest form, the model features a dividend rate that is mean-reverting around a constant fraction of the stock price. The advantage of directly specifying dynamics for the dividend rate, as opposed to the more common approach of modeling the dividend yield, is that it is easier to keep the distribution of cumulative dividends tractable. The model is non-affine but does belong to the more general class of polynomial processes, which allows us to compute all conditional moments of the stock price and the cumulative dividends explicitly. In particular, we have closed-form expressions for the prices of stock and dividend futures. Prices of stock and dividend options are accurately approximated using a moment matching technique based on the principle of maximal entropy.

Tables

Table 1. Table 1 : The second column shows market data as of 21/12/2015. DF k 𝑘 k represents the dividend futures contract with expiry on the third Friday of December ( 2015 + k ) 2015 𝑘 (2015+k) . IV stock is the Black-Scholes implied volatility of the stock option. IV dividend is the Black implied volatility of the dividend option. All data comes from Bloomberg. The last three columns show the absolute errors of the calibrated models.

		Absolute errors
	Data	$a=0.1$	$a=0.2$	$a=0.3$
DF1	115.3	0.183	0.183	0.183
DF2	108.7	0.492	0.492	0.492
DF3	105.5	1.452	1.452	1.451
DF4	100.1	0.344	0.344	0.344
DF5	95.7	0.399	0.399	0.399
DF6	92.0	0.918	0.918	0.918
DF7	89.6	0.497	0.497	0.497
DF8	87.2	0.350	0.349	0.349
DF9	84.8	0.414	0.413	0.413
DF10	84.6	1.558	1.558	1.558
IV stock	0.2295	4.095e-07	9.381e-07	9.87e-08
IV dividend	0.0491	2.001e-07	6.854e-07	9.395e-07

Table 2. Table 2 : Calibrated model parameters for three different values of a 𝑎 a .

$a$	$b$	$\beta$	$\sigma$	${\nu }_1$	$D_0$
0.1	0.0103	-0.3440	0.3621	0.0220	0.0371
0.2	0.0103	-0.3439	0.2813	0.0194	0.0371
0.3	0.0103	-0.3439	0.2614	0.0187	0.0371

Equations

$$D_t=1^{\top }{Y}_t,$$

$$\mathrm{d}{X}_t=(r{X}_t-D_t)\mathrm{d}t+\sigma (X_t-\frac{D_t}{a})\mathrm{d}{W}_t,$$

$$\mathrm{d}{Y}_t=(b{X}_t+\beta Y_t)\mathrm{d}t+\sqrt{X_t-\frac{D_t}{a}}{\left[{\nu }_1\sqrt{Y_{1,t}}\mathrm{d}{B}_{1,t},\dots ,{\nu }_d\sqrt{Y_{d,t}}\mathrm{d}{B}_{d,t}\right]}^{\top },$$

$$E=\{(x,y)\in {\mathbb{R}}^{1+d}:x>0,y\ge 0,1^{\top }y\le ax\}.$$

$$b_k+a\underset{\begin{array}{c}l=1,\dots ,d\\ l\ne k\end{array}}{\min }{\beta }_{k,l}^{-}\ge 0,\text{for all k∈{1,…,d}},$$

$$r-a-\underset{k=1,\dots ,d}{\max }(1^{\top }\beta {)}_k-\frac{1^{\top }b}{a}\ge 0.$$

$$b_k+\underset{\begin{array}{c}l=1,\dots ,d\\ l\ne k\end{array}}{\min }{\left(a{\beta }_{k,l}+\frac{{\nu }_k^2}{2}\right)}^{-}>\frac{{\nu }_k^2}{2};$$

$$r-a-\underset{k=1,\dots ,d}{\max }(\frac{{\nu }_k^2}{2a}+(1^{\top }\beta {)}_k)-\frac{1^{\top }b}{a}>0.$$

$$\mathrm{d}{\hat{Y}}_t=(\hat{b}{X}_t+\hat{\beta }{\hat{Y}}_t)\mathrm{d}t+\sqrt{X_t-\frac{D_t}{a}}{\left[{\hat{\nu }}_1\sqrt{{\hat{Y}}_{1,t}}\mathrm{d}{B}_{1,t},\dots ,{\hat{\nu }}_d\sqrt{{\hat{Y}}_{d,t}}\mathrm{d}{B}_{d,t}\right]}^{\top },$$

$$0\le {\delta }_t\le a.$$

$$\mathrm{d}{x}_t=\left(r-{\delta }_t-\frac{{\sigma }^2}{2}(1-\frac{{\delta }_t}{a}{)}^2)\right)\mathrm{d}t+\sigma (1-\frac{{\delta }_t}{a})\mathrm{d}{W}_t.$$

$${\mathbb{E}}_t\left[{\int }_t^{\infty }{\mathrm{e}}^{-r(s-t)}{D}_s\mathrm{d}s\right]=X_t.$$

$$X_t={\mathbb{E}}_t[{\mathrm{e}}^{-r(T-t)}{X}_T]+{\mathbb{E}}_t\left[{\int }_t^T{\mathrm{e}}^{-r(s-t)}{D}_s\mathrm{d}s\right],T\ge t.$$

$$\mathrm{d}{X}_t=(r{X}_t-D_t)\mathrm{d}t+\sigma (X_t-\frac{D_t}{a})\mathrm{d}{W}_t+\sqrt{X_t-\frac{D_t}{a}}\sum _{k=1}^d{\eta }_k\sqrt{Y_{k,t}}\mathrm{d}{B}_{k,t},$$

$$\mathrm{d}{X}_t=(r{X}_t-D_t)\mathrm{d}t+\sigma (X_t-\frac{D_t}{a})\mathrm{d}{W}_t,$$

$$\mathrm{d}{D}_t=(b{X}_t+\beta D_t)\mathrm{d}t+{\nu }_1\sqrt{D_t(X_t-\frac{D_t}{a})}\mathrm{d}{B}_t.$$

$$0\le b\le a(r-a-\beta ).$$

$$\frac{{\nu }_1^2}{2}<b<a(r-a-\beta )-\frac{{\nu }_1^2}{2}.$$

$$\mathrm{d}{\delta }_t=\left(b+(\beta -r){\delta }_t+{\delta }_t^2+{\sigma }^2{\delta }_t(1-\frac{{\delta }_t}{a}{)}^2\right)\mathrm{d}t-\sigma {\delta }_t(1-\frac{{\delta }_t}{a})\mathrm{d}{W}_t+{\nu }_1\sqrt{{\delta }_t(1-\frac{{\delta }_t}{a})}\mathrm{d}{B}_t.$$

$${\mathbb{E}}_t\left[{\int }_t^{\infty }{\mathrm{e}}^{-r(s-t)}{D}_s\mathrm{d}s\right]={\int }_t^{\infty }{\mathrm{e}}^{(\beta -r)(s-t)}\mathrm{d}s{D}_t=\frac{D_t}{r-\beta }.$$

$${\mathbb{E}}_t\left[{\int }_t^{\infty }{\mathrm{e}}^{-r(s-t)}{D}_s\mathrm{d}s\right]\le \frac{a}{r-\beta }{X}_t\le X_t.$$

$$C_t={\int }_0^t{D}_s\mathrm{d}s,t\ge 0,$$

$$\mathcal{G}f=(1^{\top }y,rx-1^{\top }y,(bx+\beta y{)}^{\top })\nabla f+\frac{1}{2}{\sigma }^2(x-\frac{1^{\top }y}{a}{)}^2{f}_{xx}+\frac{1}{2}\sum _{k=1}^d{\nu }_k^2{y}_k(x-\frac{1^{\top }y}{a})f_{y_k{y}_k},$$

$$\mathcal{G}{H}_n(c,x,y)=G_n{H}_n(c,x,y).$$

$${\mathbb{E}}_t[H_n(C_T,X_T,Y_T)]={\mathrm{e}}^{G_n(T-t)}{H}_n(C_t,X_t,Y_t),\forall T\ge t.$$

$$G_1=\left(\begin{array}{ccc}0& 0& 1^{\top }\\ 0& r& -1^{\top }\\ 0& b& \beta \end{array}\right).$$

$${\mathbb{E}}_t[X_T]={\mathbf{e}}_2^{\top }{\mathrm{e}}^{G_1(T-t)}{H}_1(C_t,X_t,Y_t),$$

$${\mathbb{E}}_t[C_{T_1}-C_{T_0}]={\mathbf{e}}_1^{\top }\left({\mathrm{e}}^{G_1(T_1-t)}-{\mathrm{e}}^{G_1(T_0-t)}\right)H_1(C_t,X_t,Y_t),$$

$${\mathbb{E}}_t[C_{T_1}-C_{T_0}]={\mathbf{e}}_1^{\top }{\mathrm{e}}^{G_1(T_1-t)}{H}_1(C_t,X_t,Y_t)-C_{T_0}.$$

$${\pi }_t={\mathbb{E}}_t[F(X_T)].$$

$$\begin{array}[]{ccc}\displaystyle\max_{f}&-\displaystyle\int_{0}^{\infty}f(x)\ln f(x)\,{\mathrm{d}}x&\\ \mathrm{s.t.}&\displaystyle\int_{0}^{\infty}x^{n}f(x)\,{\mathrm{d}}x=M_{n},&\quad n=0,\ldots,N,\end{array}$$

$$f^{(N)}(x)=\exp \left(-\sum _{n=0}^N{\lambda }_n{x}^n\right),$$

$${\pi }^{(N)}={\int }_0^{\infty }F(x)f^{(N)}(x)\mathrm{d}x.$$

Full text

1 Introduction

In recent years there has been an increased interest in trading dividend derivatives, in particular dividend futures. Since there also exists an active market for derivatives referencing the price of the stock paying the dividends, there is a need for derivative pricing models that can jointly price derivatives on the stock and on the dividends. Since dividend derivatives typically reference the nominal amount of dividends paid over a window of time, it seems natural to directly specify tractable dynamics for the dividend payments, or the dividend rate if dividends are paid out continuously, under a risk-neutral measure. The challenging part of this approach is to keep the stock price positive. Indeed, in absence of arbitrage and frictions such as taxes, the stock price must decrease by exactly the amount paid out as dividend, which can push the stock price in negative territory if no connection is made between the dividend and stock price dynamics. An easy solution to this problem is to model dividend yields, i.e., the fraction of the stock price that is paid out as a dividend, instead of dividends themselves. However, such a choice complicates the valuation of dividend derivatives, since their payoff now involves the product between the stock price and the dividend yield.

In this paper, we consider a stock that pays out dividends continuously222In reality, dividends are paid out discretely instead of continuously. Most liquid dividend derivatives, however, reference dividends paid by a stock index, in which case continuously paid dividends are considered an acceptable approximation. at a rate that is stochastically varying over time. The dividend rate is defined as a linear function of a multivariate factor process that belongs to the class of of polynomial diffusions, see e.g. Filipović and Larsson (2016). The drift and martingale part of the factor process and the martingale part of the stock price process are specified such that the dividend rate is non-negative and upper bounded by a constant fraction of the stock price. As a consequence, the dividend rate will go to zero as the stock price goes to zero, which guarantees the non-negativity of the stock price. The dividend yield process therefore has a zero lower bound and an upper bound that can be chosen arbitrarily large. We show that the zero boundary of the stock price is always unattainable and we provide parameter restrictions for boundary non-attainment of the factor process. The dynamics of the model becomes most intuitive with a one dimensional factor process. In this case, the dividend rate is mean-reverting around a constant fraction of the stock price. Long-term dividend futures therefore behave proportional to the stock price, while short term dividend futures have dynamics of their own. In particular, this allows to have low volatility in short-term dividend futures, while long-term dividend futures inherit the relatively high volatility of the stock price. As highlighted by Buehler (2018), this is a desirable feature for stochastic dividend models. Since the model belongs to the class of polynomial processes, we can compute all moments of the stock price and the (integrated) dividend rate in closed-form. In particular, this means we have closed-form prices for stock and dividend futures. Prices of stock and dividend options are efficiently approximated using a moment-matching technique based on the principle of maximal entropy, similarly as in Filipović and Willems (2018). In a numerical study we show that option prices are approximated accurately with a small number of moments.

We show that our model does not contain a bubble in the sense that the discounted gains process is a true martingale and that the stock price is equal to the present value of all future dividends. It is important to note that the latter is not necessarily true in an arbitrage free model. Indeed, even if the discounted gains process is a true martingale, absence of arbitrage only guarantees that the present value of all future dividends is lower than or equal to the stock price. The difference between the two is equal to the present value of the stock price at an infinite time horizon, which we show to be zero in our model. This property distinguishes our model from the one of Buehler (2018), where the dividend rate is mean-reverting around the present value of the stock price at infinity.

The literature on dividend derivative pricing is relatively scarce. Buehler et al. (2010) proposes a model where the stock price jumps at known dividend payment dates and follows log-normal dynamics in between the payment dates. The jump amplitudes are driven by an Ornstein-Uhlenbeck process such that the stock price remains log-normally distributed and the model has closed-form prices for European call options on the stock. Dividend derivatives, however, have to be priced with Monte-Carlo simulations. In order to reconcile the high volatility in the stock price with the low volatility in dividends, they use a very negative ($-95\%$) correlation between the process driving the jump amplitudes and the stock price. An important drawback of their approach is that dividend payments are not guaranteed to be non-negative. Tunaru (2018) uses a similar setup as Buehler et al. (2010), but uses a beta distribution for the jump amplitudes. The choice for a compactly supported jump distribution guarantees non-negative dividend payments. However, the diffusive noise of the stock has to be assumed independent of the jump amplitudes in order to have tractable expressions for dividend futures prices. Smoothing the dividends through a negative correlation between stock price and the jump amplitudes, as in Buehler et al. (2010), is therefore not possible. In a second approach, Tunaru (2018) directly models the cumulative dividends with a logistic diffusion model. The latter has however no guarantee to be monotonically increasing, meaning that negative dividends occur frequently. Moreover, the model must be reset on an annual basis. Option pricing is done using Monte-Carlo simulation for both methods in Tunaru (2018). Buehler (2018) decomposes the stock price as the sum of a fundamental component, representing the present value of all future dividends, and a residual bubble component. The dividends are driven by a process that mean-reverts around the bubble component. The aim of this setup is to capture the stylized fact that long-term dividend futures tend to move together with the stock price, while short-term dividend futures are much less volatile. Our model, in its simplest form, shares some similarities with the approach of Buehler (2018). However, instead of modeling the dividends as mean-reverting around a bubble component, we choose to make them mean-revert around the stock price itself, which seems more intuitive and leads directly to the desired positive correlation between long-term dividend futures and the stock price. Guennoun and Henry-Labordere (2017) consider a stochastic local volatility model for the pricing of stock and dividend derivatives. Their model guarantees a perfect fit to observed option prices, however all pricing is based on Monte-Carlo simulations. Filipović and Willems (2018) introduce a framework based on polynomial jump-diffusions to jointly price interest rate, dividend, and stock derivatives. Our model is a special case of their general framework, but is different from the Linear Jump-Diffusion model (LJD) that was used in the numerical study of Filipović and Willems (2018).

The model proposed in this paper also shares some similarities with the linear hypercube model model of Ackerer and Filipović (2019) in the context of credit risk. Specifically, they specify a survival process whose drift is a linear function of a diffusive factor process with linear drift. In order for the survival process to be positive and non-increasing, they specify the factor process such that its components are all non-negative and upper bounded by the survival process. In our setup, the dividend rate is a linear function of a diffusive factor process with linear drift, which has to be specified such that the stock price is positive and the dividend rate non-negative. The stock price, whose drift is linear in the dividend rate, therefore plays a similar role as the survival process in Ackerer and Filipović (2019), but with the important difference that the stock price has a martingale part while the survival process is absolutely continuous. This martingale part requires special care and, in particular, rules out the factor process specification of the linear hypercube model of Ackerer and Filipović (2019).

The remainder of this paper is structured as follows. Section 2 specifies the model dynamics. Section 3 discusses the pricing of stock and dividend derivatives. In Section 4 we calibrate a parsimonious model specification to market data. Section 5 presents an extension of the model with jumps in the stock price. Section 6 concludes. All proofs are collected in the appendix.

2 Model specification

Let $X_{t}$ denote the stock price process and $D_{t}$ the instantaneous dividend rate. Suppose for simplicity that interest rates are constant. Consider the following dynamics for $(X_{t},D_{t})$ under a risk-neutral measure $\mathbb{Q}$

[TABLE]

where $r\in{\mathbb{R}}$ is the short-rate, $\sigma,a>0$, $Y_{t}$ is a $d$-dimensional factor process, $d\geq 1$, $b\in{\mathbb{R}}^{d}$, $\beta\in{\mathbb{R}}^{d\times d}$, $\nu_{1},\ldots,\nu_{d}\geq 0$, and $(W_{t},B_{1,t},\ldots,B_{d,t})$ is a $(1+d)$-dimensional standard Brownian motion. The following proposition provides parameter conditions such that (2)-(3) admits a unique solution taking values in

[TABLE]

Proposition 2.1**.**

Denote by $x^{-}=\min(x,0)$. Suppose that

[TABLE]

Then for every initial value $(X_{0},Y_{0})\in E$ there exists a unique in law $E$-valued solution $(X_{t},Y_{t})$ to (2)-(3). The solution satisfies

• •

$Y_{k,t}>0$* for all $t\geq 0$ if $Y_{k,0}>0$ and*

[TABLE]

• •

$aX_{t}>\mathbf{1}^{\top}Y_{t}$* for all $t\geq 0$ if $aX_{0}>\mathbf{1}^{\top}Y_{0}$ and*

[TABLE]

We henceforth assume that the inequalities in (4) and (5) are satisfied and $(X_{0},Y_{0})\in E$. The above proposition shows in particular that we have $X_{t}>0$ and $D_{t}\geq 0$ for all $t\geq 0$. The condition in (6) can be used to bound $D_{t}$ strictly away from zero, although this is not required from an economic point of view. The diffusive term of $X_{t}$ is specified such that it vanishes at the boundary $D_{t}=aX_{t}$, which is necessary to keep the process inside $E$. The condition in (7) can be used to bound the stock price volatility strictly away from zero.

Remark 2.2**.**

Remark that the more general specification $D_{t}=\gamma^{\top}Y_{t}$, for some $\gamma\in(0,\infty)^{d}$, is equivalent to the one in (1). Indeed, if we define $\hat{Y}_{t}=CY_{t}$ with $C=\mathrm{diag}(\gamma_{1},\ldots,\gamma_{d})$, then we can write $D_{t}=\gamma^{\top}Y_{t}=\mathbf{1}^{\top}\hat{Y}_{t}$. The dynamics of $\hat{Y}_{t}$ are of the same form as the dynamics of $Y_{t}$

[TABLE]

with $\hat{b}=Cb$, $\hat{\beta}=C\beta C^{-1}$, and $\hat{\nu}_{k}=\sqrt{\gamma_{k}}\nu_{k}$, $k=1,\ldots,d$.

If we define the dividend yield $\delta_{t}=\frac{D_{t}}{X_{t}}$, then we have

[TABLE]

The dividend yield is therefore bounded from above by a parameter $a>0$ of our choice. The dynamics of the log-price $x_{t}=\log(X_{t})$ is given by

[TABLE]

The volatility of the log-price process therefore depends on the dividend yield. There is empirical evidence that dividend payout policies affect stock price volatility, see e.g. Baskin (1989), which is consistent with the dynamics of this model. In case this is not desirable, the influence of the dividend yield on the log-price volatility can be made arbitrarily small by choosing $a$ large enough.

The following proposition shows that our model does not contain a bubble in the stock price dynamics.

Proposition 2.3**.**

The discounted gains process $G_{t}={\rm e}^{-rt}X_{t}+\int_{0}^{t}{\rm e}^{-rs}D_{s}\,{\mathrm{d}}s$ is a martingale. If $\mathbf{1}^{\top}b>0$, then we have for all $t\geq 0$

[TABLE]

Equation (8) shows that the stock price is equal to the present value of all future dividends in our model. It is important to note that this is not a trivial relationship. Indeed, from no-arbitrage principles, it only follows that the present value of future dividends must be lower than or equal to the stock price, see e.g. Filipović and Willems (2018). In general, if the discounted gains process is a martingale, then

[TABLE]

A positive difference between the stock price and the present value of future dividends can be interpreted as the present value of a terminal payment at an infinite time horizon, which is difficult to reconcile with standard economic theory. Proposition 2.3 shows that, in our model, we have $\displaystyle\lim_{T\to\infty}{\mathbb{E}}_{t}[{\rm e}^{-r(T-t)}X_{T}]=0$ if $\mathbf{1}^{\top}b>0$. The derivation of this result relies on the linear drift structure of $(X_{t},Y_{t})$ and the geometry of $E$, which are a key ingredients of our model. In Example 2.5 in the next section, we illustrate a parameterization where the assumption $\mathbf{1}^{\top}b>0$ is violated and (8) does not hold.

Remark 2.4**.**

The processes $X_{t}$ and $Y_{t}$ have zero quadratic covariation. Note that this does not mean that dividends are independent of the stock price, since $X_{t}$ still enters in the drift and diffusion function of $Y_{t}$. The dynamics of $X_{t}$ can be generalized to allow for non-zero quadratic covariation with $Y_{t}$ as follows

[TABLE]

for some parameters $\eta_{k}\in{\mathbb{R}}$, $k=1,\ldots,d$. All the results in the paper are easily adjusted to accommodate this generalization.

2.1 Single factor model

For $d=1$, we obtain the following model dynamics

[TABLE]

If $\beta<0$, then $D_{t}$ is mean-reverting around $-\frac{b}{\beta}X_{t}$, with an upper bound of $aX_{t}$. The inward pointing drift conditions (4) and (5) become

[TABLE]

Boundary non-attainment is satisfied if $0<D_{0}<aX_{0}$ and

[TABLE]

The dividend yield becomes an autonomous diffusion with the following dynamics

[TABLE]

Remark that the dividend yield process is not a polynomial diffusion, due to the terms $\delta^{2}$ and $\delta^{3}$ in the drift, and $\delta^{3}$ and $\delta^{4}$ in the diffusion function. However, since $\delta_{t}$ is typically in the order of percentage points, higher powers of $\delta$ contribute relatively little to the dynamics. In particular, the dividend yield $\delta_{t}$ has approximately a linear drift $b+(\beta-r+\sigma^{2})\delta_{t}$.

We end this section with an example where the assumption in Proposition 2.3 is violated.

Example 2.5**.**

If $b=0$, then we have ${\mathbb{E}}_{t}[D_{T}]={\rm e}^{\beta(T-t)}D_{t}$ for all $T\geq t$, and (11) becomes $r-\beta\geq a>0$. The present value of future dividends is

[TABLE]

Using $aX_{t}\geq D_{t}$ and $r-\beta\geq a$ we obtain

[TABLE]

The present value of future dividends is therefore lower than or equal to the stock price, as required by absence of arbitrage. If $aX_{t}>D_{t}$ or $r-\beta>a$, then we have an example where the present value of future dividends is strictly below the stock price, i.e., ${\mathbb{E}}_{t}[{\rm e}^{-r(T-t)}X_{T}]$ does not go to zero as $T\to\infty$.

The above example shows that the presence of $X_{t}$ in the drift of $D_{t}$ is crucial for the stock price to be equal to the present value of future dividends.

3 Derivative pricing

In this section we show how to compute prices of derivatives referencing the stock price and/or the dividends paid over some time interval.

3.1 Moments of stock price and cumulative dividends

Define the cumulative dividend process as

[TABLE]

which represents all the dividends paid out over a time interval $[0,t]$. In contrast to the instantaneous dividend rate $D_{t}$, the cumulative dividend $C_{t}$ is observable in practice. The process $(C_{t},X_{t},Y_{t})$ is jointly a polynomial diffusion, so we are able to compute all conditional moments in closed form, see e.g. Filipović and Larsson (2016) for details. Let $\mathrm{Pol}_{n}$ denote the set of polynomials $p\colon{\mathbb{R}}^{2+d}\to{\mathbb{R}}$ with $1\leq\deg(p)\leq n$. Applying the infinitesimal generator ${\mathcal{G}}$ of $(C_{t},X_{t},Y_{t})$ to a twice differentiable function $f(c,x,y)$ gives

[TABLE]

where the subscripts of $f$ denote partial derivatives, $\nabla f$ the gradient of $f$, and we have omitted the function arguments for brevity. It is easily verified that ${\mathcal{G}}\mathrm{Pol}_{n}\subseteq\mathrm{Pol}_{n}$ for any $n\in{\mathbb{N}}$. Therefore, if we fix a vector of polynomial basis functions $H_{n}=(h_{1},\ldots,h_{N_{n}})^{\top}$ for $\mathrm{Pol}_{n}$, with $N_{n}=\dim(\mathrm{Pol}_{n})$, then we can find a unique matrix $G_{n}$ such that for all $(c,x,y^{\top})^{\top}\in{\mathbb{R}}^{2+d}$

[TABLE]

By definition of the infinitesimal generator, we obtain the following moment formula

[TABLE]

In particular, we can compute all the ${\mathcal{F}}_{t}$-conditional mixed moments of $(C_{T},X_{T})$ in closed form for all $T\geq t$.

Remark 3.1**.**

If one is only interested in the moments of $X_{T}$, then there is no need to augment the state with $C_{t}$, since $(X_{t},Y_{t})$ is already a polynomial diffusion on its own.

3.2 Linear derivatives

For $n=1$, we can without loss of generality choose the basis $H_{1}(c,x,y)=(c,x,y^{\top})^{\top}$. The matrix $G_{1}$ then becomes

[TABLE]

The most actively traded linear derivatives are stock futures and dividend futures. Stock futures settle on the stock price at some terminal date $T$ and dividend futures settle on the dividends paid in a time interval $[T_{0},T_{1}]$. The moment formula (12) can be used to compute prices of stock futures and dividend futures. Indeed, futures prices are given by the risk-neutral expectation of the terminal settlement price because of continuous marking-to-market, so we get

[TABLE]

where $\mathbf{e}_{k}$ denotes the $k$-th canonical basis vector in ${\mathbb{R}}^{2+d}$, $T\geq t$, and $T_{1}\geq T_{0}\geq t$. In case the reference period of the dividend futures has already started, i.e., $0\leq T_{0}\leq t\leq T_{1}$, we get

[TABLE]

Without loss of generality we can assume that $T_{0}=0$, in which case $C_{T_{0}}=0$ and $C_{t}$ is the amount of dividends already paid.

Remark that the volatility parameters $\sigma$ and $\nu_{1},\ldots,\nu_{d}$ do not enter into the prices of dividend futures, which is a consequence of the linear drift structure of $Y_{t}$. This allows us, for example, to calibrate $b$ and $\beta$ to dividend futures first, and subsequently use $\sigma$ and $\nu_{1},\ldots,\nu_{d}$ to calibrate non-linear derivatives such as stock and dividend options. The parameter $a$ also does not appear in the prices of dividend futures, however it should be noted that the value of $a$ affects the values that $b$ and $\beta$ are allowed to take, because of the inequalities (4) and (5) that we assume to be true.

3.3 Non-linear derivatives

Consider a derivative on the stock price with discounted payoff at time $T$ given by $F(X_{T})$, for some function $F$. In absence of arbitrage, its price at time $t\leq T$ is given by

[TABLE]

The probability density function of $X_{T}$ is not known explicitly, so we cannot compute $\pi_{t}$ by direct integration in general. We do however know all the moments of $X_{T}$ through the moment formula (12). In particular, if $F$ is a polynomial, then we can compute $\pi_{t}$ explicitly. If $F$ is not a polynomial, we approximate $\pi_{t}$ using the available stock price moments and the principle of maximum entropy, similarly as in Filipović and Willems (2018). Specifically, denote by $M_{n}={\mathbb{E}}_{t}[X_{T}^{n}]$, $n=1,\ldots,N$, the first $N\geq 1$ moments of the stock price. We now look for a probability density function $f$ which has the same first $N$ moments as $X_{T}$ and has maximal entropy:

[TABLE]

where we set $M_{0}=1$ so that the density integrates to one. Jaynes (1982) motivates such a choice by noting that maximizing entropy incorporates the least amount of prior information in the distribution, other than the imposed moment constraints. In this sense it is maximally noncommittal with respect to unknown information about the distribution. Straightforward functional variation with respect to $f$ gives the following unique solution $f^{(N)}$ to the optimization problem in (16)

[TABLE]

where the Lagrange multipliers $\lambda_{0},\ldots,\lambda_{N}$ have to be solved numerically from the moment constraints. Finally, we approximate $\pi_{t}$ by numerically computing the integral

[TABLE]

We can use exactly the same approach to price dividend derivatives with discounted payoff at time $T_{1}$ given by $P(C_{T_{1}}-C_{T_{0}})$, for some function $P$. All we need are the moments of $C_{T_{1}}-C_{T_{0}}$, which can be computed explicitly using the law of iterated expectations and the moment formula (12) as follows

[TABLE]

where $v_{k}$ and $w_{k}$ are the unique vectors satisfying $w_{k}^{\top}H_{k}(c,x,y)=c^{k}$ and $v_{k}^{\top}H_{n}(c,x,y)=(-c)^{n-k}w^{\top}_{k}{\rm e}^{G_{k}(T_{1}-T_{0})}H_{k}(c,x,y)$.

4 Numerical study

As an example, we calibrate the single factor model (9)–(10) using a snapshot of real market data on 21/12/2015. The stock in the calibration exercise is the Euro Stoxx 50 index, the leading blue-chip stock index in the Eurozone. The Euro Stoxx 50 index is well suited for calibrating our model since it has a liquid dividend derivatives market associated with it. The Euro Stoxx 50 dividend futures contracts are exchange traded on Eurex and reference the sum of the declared ordinary gross cash dividends (or cash-equivalent) on index constituents that go ex-dividend during a given calendar year, divided by the index divisor on the ex-dividend day. There are always ten adjacent annual contracts available for trading, with maturities every third Friday of December. We use all ten contracts in the calibration. Euro Stoxx 50 dividend options are also exchange traded on Eurex. They are options on the futures contracts, where the maturity of the option coincides with the maturity of the futures contract, which makes them effectively options on the dividends realized in a calendar year. In the calibration, we use the Black implied volatility of the option on the first dividend futures contract with at-the-money strike (i.e., strike equal to the dividend futures price). We also use the Black-Scholes implied volatility of the option on the Euro Stoxx 50 index level with maturity in three months and at-the-money strike. The prices of the dividend futures and the implied volatility of the index and dividend option are shown in the second column of Table 1. Remark that the implied volatility of the dividend option is substantially lower ($\approx 5\%$) than the implied volatility of the index option ($\approx 23\%$). This immediately shows that models with a constant dividend yield are not appropriate to price dividend derivatives, since they produce dividend payments that are roughly as volatile as the stock price itself.

We fix $r=0.01$ and $a\in\{0.1,0.2,0.3\}$. By fixing $a$, the parameter constraint in (11) becomes a linear inequality in the free parameters $b$ and $\beta$, which most optimization routines can easily deal with. In our model, $a$ determines the upper bound for the dividend yield process $\delta_{t}$. In Figure 2 we plot a proxy of the (unobservable) dividend yield $\delta_{t}$ over time, which we calculate by dividing the price of the first to expire dividend futures contract (which has a time to maturity varying between 1 day and 1 year) by the index level. We observe that between 2010 and 2016, the dividend yield proxy moves roughly between 3% and 6%, well below the three values that we consider for $a$.

We use $N=6$ moments to compute the dividend and stock option prices using the maximal entropy method described in Section 3.3. We use the gradient-free Nelder-Mead simplex optimization algorithm to find the optimal parameters $b,\beta,\sigma,\nu_{1}$, and $D_{0}$. The calibrated parameters are shown in Table 2 and the absolute pricing errors are shown in the last three columns of Table 1. The calibrated values of $b$, $\beta$, and $D_{0}$ are almost identical for different values of $a$. This is not surprising, since these parameters mainly control the term structure of dividend futures prices, and $a$ does not enter in the pricing formula (14) for the dividend futures.333Indirectly, the dividend futures prices are to some extent affected by $a$ through the inequality (11) that has to be satisfied. However, from (9) and (10) it is clear that $a$ has an impact on the volatility of the stock price and the dividend rate. Indeed, if $a$ increases, all else being equal, then the volatility of the stock price and the dividend rate increases. To offset this effect, the calibrated parameters of $\sigma$ and $\nu_{1}$ are smaller for larger $a$. From the absolute errors in Table 1, we can see that the choice of $a$ does not really matter for the quality of the calibration, since the absolute pricing errors are almost identical. The maximal relative error in the dividend futures contracts is less than $2\%$, which is a remarkably good fit for a single factor model. Figure 2 visualizes the good fit of the calibrated model with the dividend futures term structure. The option prices are matched perfectly. This is a consequence of the fact that the dividend futures prices do not depend on the martingale part of $X_{t}$ and $D_{t}$. The parameters $\sigma$ and $\nu_{1}$ therefore remain free to calibrate to the dividend and stock option.

Figure 3 plots a simulation of the dividend yield process $\delta_{t}$ over a ten years horizon with daily discretization. We use the model parameters from Table 2 with $a=0.2$, however the plot looks identical when using the calibrated parameters with $a=0.1$ or $a=0.3$. The process is roughly mean-reverting around $b/(r-\beta-\sigma^{2})=3.61\%$, which is the mean-reversion level of $\delta_{t}$ when ignoring the higher order terms in the drift. Remark that the range of values that $\delta_{t}$ takes in the simulation is very similar to the range of values in Figure 1.

Figure 4 shows the option price approximation as a function of the number of moments $N$. As a benchmark, we run a Monte-Carlo simulation with daily time steps and $10^{5}$ sample paths. In order to reduce the variance of the Monte-Carlo estimator, we use a degree one polynomial in the underlying as a control variate, where we determine the coefficients of the polynomial through a least squares regression of the simulated payoff paths on the simulated polynomial. The solid line shows the Monte-Carlo estimator and the dashed lines are the corresponding $95\%$ confidence intervals. For the stock option, the maximal entropy approximation with four moments is already within the confidence bands and the one with six moments is exactly equal to the Monte-Carlo estimator. For the dividend option, using only two moments already provides a very accurate option price approximation. The dividend option price is easier to approximate because the volatility of the dividend rate is much lower than the volatility of the stock price.

5 Extending the model with jumps

We can enrich the model dynamics by adding jumps to $X_{t}$ as follows

[TABLE]

where $D_{t}$ is the same as before and $J_{t}$ is a compensated compound Poisson process with arrival intensity $\lambda\geq 0$ and with a jump distribution $F({\mathrm{d}}z)$ that is assumed to have moments in closed-form of all orders and a support ${\mathcal{S}}\subseteq(-1,\infty)$. Remark that $D_{t}$ is still a continuous process, so that $D_{t-}=D_{t}$. Let $\tau$ denote a jump time of $J_{t}$ and suppose that $(X_{\tau-},Y_{\tau-})\in E$. From the assumption on the support of $F$, we have

[TABLE]

where equality only holds if $aX_{\tau-}=D_{\tau-}$, in which case $\Delta X_{\tau}=0$. Therefore, the results in Proposition 2.1 remain valid since $X_{t}$ behaves as in (2) in between jump times, we have $(X_{\tau},Y_{\tau})\in E$ so the process cannot jump outside of $E$, and $aX_{\tau}>D_{\tau}$ if $aX_{\tau-}>D_{\tau-}$ so jumps to the boundary are not possible.

If we denote by ${\mathcal{G}}_{J}$ the the infinitesimal generator of $(C_{t},X_{t},Y_{t})$ in the case with jumps, then we get for a twice differential function $f$

[TABLE]

where we assume that $f$ is such that the integral is finite and we have again omitted the function arguments for brevity, except in the first term of the integrand. Since the amplitudes of the jumps in $X_{t}$ depend linearly on $X_{t}$ and $Y_{t}$, it follows immediately that ${\mathcal{G}}_{J}\mathrm{Pol_{n}}\subseteq\mathrm{Pol_{n}}$. Therefore, $(C_{t},X_{t},Y_{t})$ belongs to the class of polynomial jump-diffusions and we can compute all conditional moments in closed form.

Since $X_{t}$ enters in the dynamics of $Y_{t}$, the jumps also indirectly impact the dynamics of $D_{t}$. The magnitude of the effect of a jump in $X_{t}$ on the drift of $D_{t}$ is determined by $b$. A stylized fact of index options and index dividend options is that both have a negative skew in implied volatilities. Choosing a distribution $F$ with a sufficiently negative mean produces a negative skew in implied volatilities for both stock and dividend options. We leave a calibration to option skews for future work.

Remark 5.1**.**

It is possible to introduce jumps in $Y_{t}$ as well, although one should be careful with simultaneous jumps where $D_{t}$ jumps up and $X_{t}$ jumps down in order to avoid jumping out of $E$. We do not consider this extension in this paper.

6 Conclusion

We have introduced a model for jointly pricing stock and dividend derivatives. The novelty of our approach lies in the fact that we directly model the dividend rate while guaranteeing a positive stock price. This is accomplished by upper bounding the dividend rate by a constant fraction of the stock price, so that the dividend rate goes to zero as the stock price approaches zero. The model belongs to the class of polynomial diffusions, which leads to closed-form prices for stock and dividend futures, and efficient approximations for stock and dividend options. We have calibrated a single factor model to data on dividend futures and at-the-money stock and dividend options. Future research includes calibrating the model to stock and dividend options with a range of strikes using the extension with jumps outlined in Section 5, as well as extending our framework to discrete dividend payments.

Appendix A Proofs

A.1 Proof of Proposition 2.1

Existence of an ${\mathbb{R}}^{1+d}$-valued solution follows from (Ikeda and Watanabe, 1981, Theorem IV.2.4), since the drift and dispersion coefficient of $(X_{t},Y_{t})$ satisfy a linear growth condition. It remains to show that a solution starting in $E$ also stays in $E$.

Denote by $\mu\colon E\to{\mathbb{R}}$ and $\Sigma\colon E\to{\mathbb{R}}^{d\times d}$ respectively the drift and dispersion function of $Y_{t}$, i.e.

[TABLE]

We need to verify that $\mu_{k}(x,y)\geq 0$ for $(x,y)\in E$ with $y_{k}=0$, so that the drift pushes $Y_{k,t}$ away from the zero boundary again. Using the fact that $0\leq y_{k}\leq ax$ for all $(x,y)\in E$, we have for all $(x,y)\in E$ with $y_{k}=0$ that

[TABLE]

The above inequality, together with $\Sigma_{k,l}(x,y)=0$, $l=1,\ldots,d$, for $(x,y)\in E$ with $y_{k}=0$, shows that $Y_{k,t}\geq 0$ for all $t\geq 0$ and all $k=1,\ldots,d$. Indeed, $Y_{k,t}$ starts in $E$ and has an inward pointing drift and vanishing diffusion at the boundary. Using the same argument for $aX_{t}-\mathbf{1}^{\top}Y_{t}$ instead of $Y_{k,t}$, it follows that $aX_{t}\geq\mathbf{1}^{\top}Y_{t}$ for all $t\geq 0$. As a consequence we also have $X_{t}\geq 0$.

In order to prove the non-attainment of the zero lower boundary of $X_{t}$, we use a stochastic comparison argument. Define the process $Z_{t}=-\log X_{t}$ if $X_{t}>0$ and $Z_{t}=\infty$ if $X_{t}=0$. Define the process $\tilde{Z}_{t}$ through the following SDE

[TABLE]

From Theorem 1.3 in Hajek (1985) we get for all $c\in{\mathbb{R}}$ and $t>0$

[TABLE]

Since $\displaystyle\lim_{c\to\infty}P(\tilde{Z}_{t}\geq c)=0$, we have $\displaystyle\lim_{c\to\infty}P(Z_{t}\geq c)=0$ and therefore $X_{t}>0$ a.s.

We use Theorem 5.7(i) in Filipović and Larsson (2016), which we restate below for completeness, to study boundary attainment of $Y_{t}$.

Theorem A.1** **(Theorem 5.7(i) Filipović and

Larsson (2016)).

Denote by ${\mathcal{G}}$ the infinitesimal generator and by $m(x,y)$ the diffusion function of $(X_{t},Y_{t})$. Let $p(x,y)$ be a polynomial and let $h(x,y)$ be a vector of polynomials such that $m(x,y)\nabla p(x,y)=h(x,y)p(x,y)$ for all $(x,y)\in{\mathbb{R}}^{1+d}$. If there exists a neighborhood $U$ of $E\cap\{p=0\}$ such that for all $(x,y)\in E\cap U$

[TABLE]

then $p(X_{t},Y_{t})>0$ for all $t>0$.

First, we derive conditions such that $Y_{k,t}>0$. For $p(x,y)=y_{k}$, we have

[TABLE]

with the non-zero element in the $(k+1)$-th component. For some $\epsilon>0$, consider the following neighborhood of $E\cap\{p=0\}$

[TABLE]

For $(x,y)\in E\cap U=\{(x,y)\in{\mathbb{R}}^{1+d}\colon x>0,y_{k}\leq\epsilon,ax\geq\mathbf{1}^{\top}y,y\geq 0\}$ we have

[TABLE]

If $2\beta_{k,k}+\nu_{k}^{2}/a\geq 0$, then (19) is non-negative if $2b_{k}-\nu_{k}^{2}+\min_{l\neq k}(2a\beta_{k,l}+\nu_{k}^{2})^{-}\geq 0$. If $2\beta_{k,k}+\nu_{k}^{2}/a<0$, then we can always find an $\epsilon>0$ such that (19) is non-negative if $2b_{k}-\nu_{k}^{2}+\min_{l\neq k}(2a\beta_{k,l}+\nu_{k}^{2})^{-}>0$.

Finally, we derive conditions such that $aX_{t}>\mathbf{1}^{\top}Y_{t}$. For $p(x,y)=ax-\mathbf{1}^{\top}y$ we have

[TABLE]

For some $\epsilon>0$, consider the following neighborhood of $E\cap\{p=0\}$

[TABLE]

For $(x,y)\in E\cap U=\{(x,y)\in{\mathbb{R}}^{1+d}\colon x>0,0\leq ax-\mathbf{1}^{\top}y\leq\epsilon,y\geq 0\}$ we have

[TABLE]

where the last line follows from $ax-\epsilon\leq\mathbf{1}^{\top}y\leq ax$. If $2a-\sigma^{2}+\max_{k=1,\ldots,d}(\nu_{k}^{2}/a+2\mathbf{1}^{\top}\beta_{k})\geq 0$, then $2{\mathcal{G}}p-h^{\top}\nabla p\geq 0$ on $E\cap U$ if

[TABLE]

If $2a-\sigma^{2}+\max_{k=1,\ldots,d}(\nu_{k}^{2}/a+2\mathbf{1}^{\top}\beta_{k})<0$, then we can always find an $\epsilon>0$ such that $2{\mathcal{G}}p-h^{\top}\nabla p\geq 0$ on $E\cap U$ if

[TABLE]

For uniqueness in law of the solution $(X_{t},Y_{t})$, note that $\frac{Y_{t}}{X_{t}}$ is an autonomous diffusion with $0\leq\frac{Y_{t}}{X_{t}}\leq a$ for all $t\geq 0$. A straightforward application of Itô’s lemma shows that the process $(\log(X_{t}),\frac{Y_{t}}{X_{t}})$ has a uniformly bounded drift and diffusion function, so that uniqueness in law for $(\log(X_{t}),\frac{Y_{t}}{X_{t}})$, and therefore for $(X_{t},Y_{t})$, follows from (Ikeda and Watanabe, 1981, Theorem IV.3.3).

A.2 Proof of Proposition 2.3

To proof that $G_{t}$ is a martingale, we can use Novikov’s condition. An application of Itô’s lemma gives

[TABLE]

which shows that $G_{t}$ is a local martingale. Since the volatility of $\log(G_{t})$ is uniformly bounded,

[TABLE]

Novikov’s condition is trivially satisfied and we conclude that $G_{t}$ is a martingale.

Remark A.2**.**

The process $G_{t}$ represents the discounted value of a trading strategy of a long position in the stock and investing all the dividends in the risk-free account. Alternatively, we could also re-invest all the dividends in the stock itself. This strategy has a discounted value process $G^{\ast}_{t}={\rm e}^{-rt+\int_{0}^{t}\delta_{s}\,{\mathrm{d}}s}X_{t}$, which is again a martingale by Novikov’s condition.

Next, we show that the present value of all future dividends is equal to the stock price. Define $\tilde{X}_{t}={\rm e}^{-rt}X_{t}$ and $\tilde{Y}_{t}={\rm e}^{-rt}Y_{t}$. The dynamics of $\tilde{X}_{t}$ and $\tilde{Y}_{t}$ becomes

[TABLE]

Taking conditional expectations and denoting $f(s)={\mathbb{E}}_{t}[\tilde{X}_{s}]$ and $g(s)={\mathbb{E}}_{t}[\tilde{Y}_{s}]$, $s\geq t$, gives the following linear first order ODE

[TABLE]

Using the properties of $E$, we have that $f(s)>0$, $g(s)\geq 0$, and $\mathbf{1}^{\top}g(s)\leq af(s)$ for all $s\geq t$. In particular, we have that $f$ is a non-increasing function and hence $f$ and $g$ are uniformly bounded

[TABLE]

Moreover, all derivatives of $f$ and $g$ are uniformly bounded as well since

[TABLE]

for all $n\in{\mathbb{N}}$. Since $f$ is a non-increasing positive function, we have $\displaystyle\lim_{s\to\infty}f(s)=\xi\in[0,f(t)]$. Since $f^{\prime\prime}$ is bounded, $f^{\prime}$ must be uniformly continuous. By Barbalat’s lemma (see e.g., Lemma 8.2 in Khalil (2002)) we therefore have that $\displaystyle\lim_{s\to\infty}f^{\prime}(s)=0$. Since $f^{\prime}(s)=-\mathbf{1}^{\top}g(s)$ and $g(s)\geq 0$, we also have $\displaystyle\lim_{s\to\infty}g(s)=0$ componentwise. Taking the limit of $f^{\prime\prime}$ gives

[TABLE]

Since $f^{\prime\prime\prime}$ is bounded, $f^{\prime\prime}$ must be uniformly continuous, and by Barbalat’s lemma we have $\displaystyle\lim_{s\to\infty}f^{\prime\prime}(s)=0$. Since $\mathbf{1}^{\top}b>0$ by assumption, we must have $\xi=0$. This concludes the proof since

[TABLE]