Find Conditional Density for Joint Continuous

Probability Distributions I

B.R. Martin , in Statistics for Physical Science, 2012

3.3.2 Marginal and Conditional Distributions

We may also be interested in the density function of a subset of variables. This is called the marginal density function $f^M$ , and in general for a subset $x_1(i=1,2,\dots ,m<n)$ of the variables is given by integrating the joint density function over all the variables other than $x_1,x_2,\dots ,x_m$ . Thus,

(3.20a) $f^M(x_1,x_2,\dots ,x_m)=\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}\cdots \underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}f(x_1,x_2,\dots ,x_m,x_{m+1},\dots ,x_n)\prod _{i=m+1}^n\mathrm{d}{x}_i.$

In the case of two variables, we may be interested in the density function of x regardless of the value of y, or the density function of y regardless of x. For example, the failure rate of a resistor may be a function of its operating temperature and the voltage across it, but in some circumstances we might be interested in just the dependence on the former. In these cases (3.20a) becomes

(3.20b) $f^M(x)=\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}f(x,y)\text{d}y\text{and}{f}^M(y)=\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}f(x,y)\mathrm{d}x.$

These density functions correspond to the normalized histograms obtained by projecting a scatter plot of x and y onto one of the axes. This is illustrated in Fig. 3.4, again using the data of Fig. 1.3(b).

FIGURE 3.4. Normalized histograms obtained by projecting the data of Fig. 1.3(b) onto the x and y axes, together with the corresponding marginal probability density functions $f^M(x)$ and $f^M(y)$ .

We can also define the multivariate conditional density function of the random variables $x_1(i=1,2,\dots ,m<n)$ by

(3.21) $f^C(x_1,x_2,\dots ,x_m|x_{m+1},x_{m+2},\dots ,x_n)\equiv \frac{f(x_1,x_2,\dots ,x_n)}{f(x_{m+1},x_{m+2},\dots ,x_n)}.$

Again, if we consider the case of two variables x and y, the probability for y to be in the interval $(y,y+\text{d}y)$ with any x (event B), given that x is in the interval $(x,x+\text{d}x)$ with any y (event A), is

$P[B|A]=\frac{P[A\cap B]}{P[A]}=\frac{f(x,y)\mathrm{d}x\mathrm{d}y}{f^M(x)\mathrm{d}x},$

where $f^M(x)$ is the marginal density function for x. The conditional density function for y given x, is thus

(3.22a) $f^C(y|x)=\frac{f(x,y)}{f^M(x)}=\frac{f(x,y)}{\int f(x,y^{\prime })\mathrm{d}{y}^{\prime }}.$

This is the density function of the single random variable y where x is treated as a constant. It corresponds to projecting the events in a band dx centered at some value $x$ onto the y axis and renormalizing the resulting density so that it is unity when integrated over y. The form of $f^C(y|x)$ will therefore vary as different values of x are chosen.

The conditional density function for x given y is obtained from (3.22a) by interchanging x and y, so that

(3.22b) $f^C(x|y)=\frac{f(x,y)}{f^M(y)}=\frac{f(x,y)}{\int f(x^{\prime },y)\mathrm{d}{x}^{\prime }},$

and combining these two equations gives

(3.22c) $f^C(x|y)=\frac{f^C(y|x)f^M(x)}{f^M(y)},$

which is Bayes' theorem for continuous variables.

We can use these definitions to generalize the law of total probability (2.6b) to the case of continuous variables. Using conditional and marginal density functions we have

(3.23) $f(x,y)=f^C(y|x)f^M(x)=f^C(x|y)f^M(y),$

so the marginal density functions may be written as

$f^M(y)=\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}{f}^C(y|x)f^M(x)\mathrm{d}x$

and

$f^M(x)=\underset{-\mathrm{\infty }}{\overset{+\mathrm{\infty }}{\int }}{f}^C(x|y)f^M(y)\text{d}y.$

With more than one random variable we also have to consider the question of statistical independence (by analogy with the work of Chapter 2). If the random variables may be split into groups such that their joint density function is expressible as a product of marginal density functions of the form

$f(x_1,x_2,\dots ,x_n)=f_1^M(x_1,x_2,\dots ,x_i)f_2^M(x_{i+1},x_{i+2},\dots ,x_k)\dots f_n^M(x_{l+1},x_{l+2},\dots ,x_n),$

then the sets of variables

$(x_1,x_2,\dots ,x_i);(x_{i+1},x_{i+2},\dots ,x_k);\dots ;(x_{l+1},x_{l+2},\dots ,x_n),$

are said to be statistically independent, or independently distributed. So two random variables x and y are independently distributed if

(3.24) $f(x,y)=f^M(x)f^M(y).$

It follows from (3.22) that in this case the conditional density function of one variate does not depend on knowledge about the other variate.

EXAMPLE 3.7

The joint mass function for two discrete variables x and y is given by

$f(x,y)=\{\begin{array}{l}k(2x+3y)0\le x\le 3,0\le y\le 2\\ 0\mathit{otherwise}\end{array},$

where k is a constant. Find: (a) the value of k, (b) $P[x\ge 2,y\le 1]$ , and (c) the marginal density of x.

The mass function may be tabulated as below.

$\begin{array}{cccccc}& & & y& & \\ x& & 0& 1& 2& \text{total}\\ 0& & 0& 3k& 6k& 9k\\ 1& & 2k& 5k& 8k& 15k\\ 2& & 4k& 7k& 10k& 21k\\ 3& & 6k& 9k& 12k& 27k\\ \text{total}& & 12k& 24k& 36k& 72k\end{array}$

(a)

The normalization condition is $\sum _{x,y}f(x,y)=1$ , so $k=1/72$ .

(b)

$\begin{array}{c}P[x\ge 2,y\le 1]=P[x=2,y=1]+P[x=2,y=0]+P[x=3,y=1]+P[x=3,y=0]\\ =7k+4k+9k+6k=26k=13/36.\end{array}$

(c)

The marginal probability of x is

$P[x]=\sum _{y}P[x,y]=\{\begin{array}{l}9k=3/24x=0\\ 15k=5/24x=1\\ 21k=7/24x=2\\ 27k=9/24x=3\end{array}$

EXAMPLE 3.8

If $f(x,y)$ is the joint density function of two continuous random variables x and y, defined by

$f(x,y)=\{\begin{array}{l}{e}^{-(x+y)}x,y\ge 0\\ 0otherwise\end{array},$

find their conditional distribution.

From (3.22b),

$f^C(x|y)=\frac{f(x,y)}{f^M(y)},$

where the marginal density of y is given from (3.20b) as

$f^M(y)=\underset{0}{\overset{\mathrm{\infty }}{\int }}f(x,y)\mathrm{d}x=e^{-y}{[-e^{-x}]}_0^{\mathrm{\infty }}=e^{-y}.$

Thus

$f^C(x|y)=\frac{e^{-(x+y)}}{e^{-y}}=e^{-x}.$

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780123877604000032

Conceptual Econometrics Using R

Kris Boudt , ... Eric Zivot , in Handbook of Statistics, 2019

6 BIP and GAS MGARCH models

The estimated covariance updating equation of MGARCH models is mostly used as a filter to predict the conditional covariance based on the observed return series. The already discussed MGARCH models tend to use the same filter irrespective of the shape of the distribution function. This typically leads to a large spike in the conditional covariance prediction following an extreme return realization. In order to dampen the effect of outliers on covariance predictions, Boudt and Croux (2010) and Boudt et al. (2013) recommend to use robust MGARCH filters that have the property of Bounded Innovation Propagation (BIP). The resulting model is called BIP-MGARCH for which they propose robust M-estimators under the assumption of elliptical innovations. The BIP-BEKK model corresponding to the BEKK model in (5) is given by:

(63) ${\mathit{H}}_t={\mathit{C}}^{\prime }\mathit{C}+\sum _{k=1}^K\sum _{j=1}^q{\mathit{A}}_{\mathit{jk}}^{\prime }{\tilde{\varepsilon }}_{t-j}{\tilde{\varepsilon }}_{t-j}^{\prime }{\mathit{A}}_{\mathit{jk}}+\sum _{k=1}^K\sum _{j=1}^p{\mathit{B}}_{\mathit{jk}}^{\prime }{\mathit{H}}_{t-j}{\mathit{B}}_{\mathit{jk}},$

with

(64) ${\tilde{\varepsilon }}_t={\varepsilon }_t\sqrt{w\left({\varepsilon }_t^{\prime }{\mathit{H}}_t^{-1}{\varepsilon }_t\right)}.$

The weight function w(⋅) must be such that the effect of ɛ _t on H _t is bounded. Boudt and Croux (2010) use the following weight function

(65) $w\left(z\right)=\left\{\begin{array}{lll}1& & \text{if}z\le c_1\\ 1-{\left(1-c_1/z\right)}^3& & \text{if}{c}_1<z\le c_2\\ \left(c_2/z\right)\left(1-{\left(1-c_1/c_2\right)}^3\right)& & \text{else}.\end{array}\right.$

They set the parameters c ₁ and c ₂ equal the 99% and 99.9% quantile of the distribution of the squared Mahalanobis Distances (MD) ɛ _t ^′ H _t ^{−   1}ɛ _t . The left panel of Fig. 1 shows this weight function in case the conditional distribution of ɛ _t is bivariate Normal (top plot) or Student t ₄ (bottom plot). Note that only the observations with an extremely large MD are downweighted and that the weighting depends on the distributional assumption. In the right panel we also plot the function w(z)z, which is of interest since if z is the squared MD of ${\tilde{\varepsilon }}_t$ , then w(z)z is the squared MD of ɛ _t . Note that the downweighting is such that the function w(z)z is nondecreasing and bounded by w(c ₂)c ₂. The smoothness of the weight function is needed to avoid numerical problems in the parameter estimation.

Fig. 1. Plot of the functions w(z) and w(z)z used in the bivariate BIP-BEKK model. The parameters c ₁ and c ₂ equal to the 99% and 99.9% quantiles of the squared MD under Gaussian (upper panel) and Student t ₄ innovations (lower panel).

The BIP-BEKK model can be seen as an ad hoc robustification of the BEKK filters. Similar BIP-DCC filters were proposed in Boudt et al. (2013). An elegant alternative to take the shape of the distribution function into account is the class of Generalized Autoregressive Score (GAS) and Dynamic Score models, proposed by Creal et al. (2013) and Harvey (2013) at Vrije Universiteit Amsterdam and Cambridge University, respectively. They developed a general framework to specifying the time-varying parameters of a conditional distribution function. The key feature of their framework is that the score of the conditional density function is used as the driver of the time variation in the parameters, making it possible to obtain the likelihood in closed form through a standard prediction error decomposition.

More formally, suppose that the variable of interest is r _t with conditional density function f(…). The conditional density depends on a vector of time-varying parameters denoted by θ _t   ∈ Θ  ⊆   ℜ ^J . It typically contains the unique elements in the conditional covariance matrix H _t , but may also consist of location and shape parameters, among others. Usually, the parameter space of θ _t is restricted by various conditions, such as the requirement of positive definiteness for the conditional covariance. The standard solution under the GAS framework is to work with parameter transformations such that the parameter of interest θ _t is in the parameter space. More precisely, when the unrestricted parameter vector is denoted by ${\tilde{\theta }}_{\mathbf{t}}\in {\mathrm{\Re }}^{\mathbf{J}}$ , the GAS model uses a link function Λ(⋅) to map the transformed parameter ${\tilde{\theta }}_{\mathbf{t}}\in {\mathrm{\Re }}^{\mathbf{J}}$ into the parameter of interest θ _t . The evolution in the time-varying parameter vector ${\tilde{\theta }}_{\mathbf{t}}$ is driven by the scaled score of the conditional density function, defined as:

${\mathbf{s}}_{\mathbf{t}}={\mathbf{S}}_{\mathbf{t}}\frac{\partial log\mathbf{f}\left({\mathbf{r}}_{\mathbf{t}}\left|{\theta }_{\mathbf{t}}\right.\right)}{\partial {\tilde{\theta }}_{\mathbf{t}}},$

where the matrix S _t is a J  × J positive definite scaling matrix known at time t. ^w The quantity s _t indicates the direction to update the vector of parameters from θ _t , to θ _{t   +   1} , acting as a steepest ascent algorithm for improving the model's local fit given the current parameter position. This updating procedure resembles the well-known Newton–Raphson algorithm.

The updating equation for the unconstrained parameter vector is given by a linear function of the score, together with an autoregressive component:

(66) ${\tilde{\theta }}_{\mathbf{t}\mathbf{+}\mathbf{1}}=\kappa +\mathbf{A}{\mathbf{s}}_{\mathbf{t}}+\mathbf{B}{\tilde{\theta }}_{\mathbf{t}},$

where κ, A, and B are matrices of coefficients with proper dimensions.

A general implementation of univariate and multivarate GAS models can be found in the $\mathsf{R}$ package GAS (Ardia et al., 2018b; Catania et al., 2017). As a simple illustration we show in code Snippet 3 how one can define a multivariate GAS model in the GAS package and estimate it. Table 4 provides the methods and functions available for working with the model which include forecasting, filtering, simulation, and inference, among others. We refer the reader to the documentation of the package for more details.

Snippet 3. GAS example.

Table 4. GAS multivariate GAS model functions and methods

Functions/methods	Description	Input classes
MultiGASSpec	Model specification	NA
MultiGASFit	Model estimation	1
getFilteredParameters	1-ahead ahead filtering	1
ConfidenceBands	Confidence bands for the filtered parameters	1
getMoments	Extract conditional moments	1, 3, 5
getForecast	Extract parameter forecast	3, 5
LogScore	Extract log scores	3, 5
MultiGASFor	1- to h-ahead forecasts	2
MultiGASSim	Simulation	2
residuals	Conditional mean equation residuals	2
coef	Coefficients of model	2
show	Summary of output	3,4,5
summary	Summary of output	2
MultiGASRoll	Rolling estimation/forecasting	1

Note: The table provides a list of the methods and functions for working with multivariate GAS models in the GAS package. The input classes are as follows: 1   =   mGASSpec, 2   =   mGASFit, 3   =   mGASFor, 4   =   mGASSim, 5   =   mGASRoll.

The approach of specifying the time variation in all the distribution parameters jointly as a function of the conditional score is unfortunately not feasible in large-scale applications due to a curse of dimensionality. Creal et al. (2011) acknowledge this shortcoming and propose to use a time-varying copula specification in order to model the variances separately from the correlations. As in the DCC model of Engle (2002), it is then straightforward to combine the conditional variances and correlations into an estimate for the conditional covariance matrix H _t . We illustrate this copula-approach next in the case of a GAS variance model assuming a Student t marginal distribution, and a GAS correlation model under the assumption of a bivariate t-copula specification.

In terms of modeling the conditional variance dynamics, we focus here on the case where a Student t distribution is assumed. Several GAS models exist under this framework. The Beta-t-EGARCH model introduced by Harvey et al. (2008) uses the exponential function as link function. The Beta-t-GARCH model of Harvey et al. (2008) and the t-GAS model of Creal et al. (2013) use no transformation. The latter then leads to a GAS volatility model that is the close to the GARCH(1,1) model and for which the estimates are published on https://vlab.stern.nyu.edu. Under this model, the conditional variance for asset i with zero mean and ν _i degrees of freedom is given by:

(67) $h_{\mathit{ii},t+1}={\omega }_i+{\alpha }_i\frac{{\nu }_i+3}{{\nu }_i}\left(\frac{{\nu }_i+1}{{\nu }_i-2+{\mathrm{ϵ}}_{i,t}^2/h_{\mathit{ii},t}}{ϵ}_{i,t}^2-h_{\mathit{ii},t}\right)+{\beta }_i{h}_{\mathit{ii},t}.$

Note that if ν _i   =   ∞   , the GAS-t volatility model collapses to a traditional GARCH model and the score has a quadratic impact on the conditional variance. This can be seen as well in Fig. 2, where we show the scaled score for various value of ν _i . Note that, the more fat-tailed the return distribution is, the more their effect on future variance is dampened due to to downweighting when ν _i is small. This is desired, since in case of a fat-tailed distribution, a large value of ϵ _i,t /h _i,t ^1/2 may as well be a tail realization and thus does not necessitate a substantial increase in the conditional variance.

Fig. 2. Scaled score used as driver in the t-GAS conditional variance model for various values of ν _i when h _ii,t   =   1.

We now discuss the GAS correlation model under the assumption of a bivariate t-copula specification. Besides the use of hyperspherical coordinate transformation, Creal et al. (2011) discuss also the approach to decompose the correlation matrix R _t as Δ _t ^{−  1} Q _t Δ _t ^{−  1} , where Q _t is a symmetric positive definite matrix and Δ _t is a diagonal matrix whose nonzero elements equal the square root of the diagonal elements of Q _t . They then use the GAS framework in (66) to obtain a score-driven calibration of the time variation in   vech   (Q _t ). The score takes the shape of the student t copula into account, together with the surprise in the realized correlations compared to the predicted ones. We refer to Creal et al. (2011) for the general specification, and limit ourselves here to summarizing the discussion in Creal et al. (2011) regarding the bivariate case with fixed unit variance, correlation parameter ρ _t , and Gaussian copula. Then, the scaled score is given by:

(68) $s_t=\frac{2}{{\left(1-{\rho }_t^2\right)}^2}\left[\left(1+{\rho }_t^2\right)\left(y_{1,t}{y}_{2,t}-{\rho }_t\right)-{\rho }_t\left(y_{1,t}^2+y_{2,t}^2-2\right)\right].$

The first term increases the conditional correlation when the realized correlation exceeds the conditional correlation, while the second term attenuates this correlation increase in case of a large dispersion in the input vector. In fact, the correlation signal of (1,   1) is much stronger than (1/4,   4), even though their cross product is the same. Boudt et al. (2012) present similar expressions for the scaled score in case of the t-copula. The most important change is that, alike in the univariate case, the input data are then also downweighted when they have a large Mahalanobis distance. The more fat-tailed the t-copula is, the larger the downweighting in order to compensate for the fact that large deviations of realized correlations from predicted ones may be tail events and thus do not necessarily imply changes in the conditional correlation.

The application of the GAS framework to create MGARCH models is still an active field of research, as can be seen on the overview website at http://www.gasmodel.com/. A large number of authors have worked out the GAS dynamics in case of a more flexible distribution function, like the univariate skewed t distribution (see, e.g., Harvey and Sucarrat (2014) and Ardia et al. (2018a)), the generalized hyperbolic skewed t distribution (Lucas et al., 2014), or the Wishard distribution (Gorgi et al., 2019) and general finite mixture of distributions (Catania, 2016). Others have generalized the GAS specification to account for regime switches (Boudt et al., 2012; Catania, 2018). In terms of attempts at obtaining GAS models that are feasible in high dimensions, we refer the reader to Boudt et al. (2012) and Lucas et al. (2017) for an analysis assuming (block) equicorrelation, and to Creal et al. (2011, 2014) for the use of dynamic factors under the GAS framework.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/S016971611930001X

Classical and Quantum Information Theory

Dan C. Marinescu , Gabriela M. Marinescu , in Classical and Quantum Information, 2012

Proof.

We assume that the discrete random variable, X, has a probability distribution, p_X(x), and call X the range of values of X. We base our estimation of X on the random variable, Y, which is related to X by the conditional density function, p_Y|X(y|x); once we know Y, we estimate X as a function of Y,

=g(Y). Clearly, X, Y,

form a Markov chain, X ↦ Y ↦

. Now we introduce a binary random variable, A, defined as

$A=\{\begin{array}{c}1if\tilde{X}\ne X,\\ 0if\tilde{X}=X,\end{array}$

We apply the chain rule for entropy and write

$H(A,X|Y)=H(X\left|Y\right)+H\left(A\right|X,Y)=H(X\left|Y\right)$

.

Indeed, H (A|X, Y) = 0 as A is completely determined once we know X and Y. We apply again the chain rule for entropy to express H (A,X|Y) differently:

$H(A,X|Y)=H(A\left|Y\right)+H\left(X\right|A,Y)$

.

A is a binary random variable, and Prob(A = 0) = 1 − p_err and Prob(A = 1) = p_err . Thus, H (A) = H(p_err ); also, H(A|Y) ≤ H(A), since conditioning reduces the entropy. Now,

$H\left(X\right|A,Y)=\text{Prob}(A=0)\cdot H(X|Y,A=0)+\text{Prob}(A=1)\cdot H\left(X\right|Y,A=1)$

.

We recognize that H(X|Y, A = 0) = 0 because X =

when A = 0; |χ| is the cardinality of the target set of X; thus, H(X|Y, A = 1) = log(|χ| − 1), and we conclude that

$H\left(X\right|A,Y)\ge H(p_{err})+(1-p_{err})\cdot 0+p_{error}\log (\left|\mathcal{X}\right|-1)=H(p_{err})+p_{error}\log (\left|\mathcal{X}\right|-1)$

.

Combining the two expressions of H(X|A, Y), we obtain the Fano inequality:

$H\left(p_{err}\right)+p_{error}\log \left(\right|\mathcal{X}|-1)\ge H\left(X\right|Y)$

.

But H(p_err ) ≤ 1 and |χ|, the cardinality of the set of sample values of the random variable X, is expected to be large, |χ| − 1 ≈ |χ|; therefore, weaker forms of Fano's inequality are

$1+p_{error}log\left(\left|\chi \right|\right)\ge H\left(X|Y\right)\text{or}{p}_{err}\ge \frac{H\left(X|Y\right)-1}{log\left(\left|\chi \right|\right)}.$

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780123838742000035

Poisson Processes

Mark A. Pinsky , Samuel Karlin , in An Introduction to Stochastic Modeling (Fourth Edition), 2011

Problems

5.4.1

Let W ₁, W ₂, … be the event times in a Poisson process {X (t); t ≥ 0} of rate λ. Suppose it is known that X (1) = n. For k ≤ n , what is the conditional density function of W ₁, …, W_{k −1}, W_{k + 1}, …, W_n , given that W_k = w?

5.4.2

Let {N (t); t ≥ 0} be a Poisson process of rate λ, representing the arrival process of customers entering a store. Each customer spends a duration in the store that is a random variable with cumulative distribution function G. The customer durations are independent of each other and of the arrival process. Let X (t) denote the number of customers remaining in the store at time t, and let Y (t) be the number of customers who have arrived and departed by time t. Determine the joint distribution of X (t) and Y (t).

5.4.3

Let W ₁, W ₂,… be the waiting times in a Poisson process {X (t); t ≥ 0} of rate λ. Under the condition that X (1) = 3, determine the joint distribution of U = W ₁/W ₂ and V = (1 − W ₃)/(1 − W ₂).

5.4.4

Let W ₁, W ₂, … be the waiting times in a Poisson process {X (t); t ≥ 0} of rate λ. Independent of the process, let Z ₁, Z ₂, … be independent and identically distributed random variables with common probability density function f (x), 0 < x < ∞. Determine Pr{Z > z}, where

$Z=min\{W_1+Z_1,W_2+Z_2,\dots \}.$

5.4.5

Let W ₁, W ₂, … be the waiting times in a Poisson process {N (t); t ≥ 0} of rate λ. Determine the limiting distribution of W ₁, under the condition that N (t) = n as n → ∞ and t → ∞ in such a way that n/t = β > 0.

5.4.6

Customers arrive at a service facility according to a Poisson process of rate λ customers/hour. Let X (t) be the number of customers that have arrived up to time t. Let W ₁, W ₂, … be the successive arrival times of the customers.

(a)

Determine the conditional mean E [W ₁| X (t) = 2].

(b)

Determine the conditional mean E [W ₃| X (t) = 5].

(c)

Determine the conditional probability density function for W 2, given that X (t) = 5.

5.4.7

Let W ₁, W ₂, … be the event times in a Poisson process {X (t); t ≥ 0} of rate λ, and let f (w) be an arbitrary function. Verify that

$E\left[\sum _{i=1}^{X\left(t\right)}f\left(W_i\right)\right]=\lambda \underset{0}{\overset{t}{\int }}f\left(w\right)\text{d}w.$

5.4.8

Electrical pulses with independent and identically distributed random amplitudes ξ₁, ξ₂, … arrive at a detector at random times W ₁, W ₂, … according to a Poisson process of rate λ. The detector output θ _k (t) for the k th pulse at time t is

${\theta }_k\left(t\right)=\{\begin{array}{cc}0& \text{for}t<W_k,\\ {\xi }_{k}exp\{-\alpha \left(t-W_k\right)\}& \text{for}t\ge W_k.\end{array}$

That is, the amplitude impressed on the detector when the pulse arrives is ξ _k , and its effect thereafter decays exponentially at rate α. Assume that the detector is additive, so that if N (t) pulses arrive during the time interval [0, t], then the output at time t is

$Z\left(t\right)=\sum _{k=1}^{N\left(t\right)}{\theta }_k\left(t\right).$

Determine the mean output E [Z (t)] assuming N (0) = 0. Assume that the amplitudes ξ₁, ξ₂, … are independent of the arrival times W ₁, W ₂,….

5.4.9

Customers arrive at a service facility according to a Poisson process of rate λ customers per hour. Let N (t) be the number of customers that have arrived up to time t, and let W ₁, W ₂, … be the successive arrival times of the customers. Determine the expected value of the product of the waiting times up to time t. (Assume that W ₁ W ₂ ··· W_{N (t)} = 1 when N (t) = 0.)

5.4.10

Compare and contrast the example immediately following Theorem 5.7, the shot noise process of Section 5.4.1, and the model of Problem 4.8. Can you formulate a general process of which these three examples are special cases?

5.4.11

Computer Challenge Let U ₀, U ₁, … be independent random variables, each uniformly distributed on the interval (0, 1). Define a stochastic process {S_n } recursively by setting

$S_0=0\text{and}{S}_{n+\text{1}}=U_n\left(1+S_n\right)\text{for}n>0.$

(This is an example of a discrete-time, continuous-state, Markov process.) When n becomes large, the distribution of S_n approaches that of a random variable S = S _∞, and S must have the same probability distribution as U (1 + S), where U and S are independent. We write this in the form

from which it is easy to determine that $E[S]=1,\text{Var}[S]=\frac{\text{1}}{\text{2}}$ , and even (the Laplace transform)

$E\left[{\text{e}}^{-\theta S}\right]=exp\left\{-\underset{0+}{\overset{\theta }{\int }}\frac{1-{\text{e}}^{-u}}{u}\text{d}u\right\},\theta >0.$

The probability density function f (s) satisfies

$\begin{array}{l}f\left(s\right)=0\text{for}s\le 0,\text{and}\\ \frac{\text{d}f}{ds}=\frac{1}{s}f\left(s-1\right),\text{for}s>0.\end{array}$

What is the 99th percentile of the distribution of S? (Note: Consider the shot noise process of Section 5.4.1. When the Poisson process has rate λ= 1 and the impulse response function is the exponential h (x) = exp {– x}, then the shot noise I (t) has, in the limit for large t, the same distribution as S.)

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780123814166000058

Conditional Probability and Conditional Expectation

Sheldon M. Ross , in Introduction to Probability Models (Twelfth Edition), 2019

3.3 The Continuous Case

If X and Y have a joint probability density function $f(x,y)$ , then the conditional probability density function of X, given that $Y=y$ , is defined for all values of y such that $f_Y(y)>0$ , by

$f_{X|Y}(x|y)=\frac{f(x,y)}{f_Y(y)}$

To motivate this definition, multiply the left side by dx and the right side by $(dxdy)/dy$ to get

$\begin{array}{cc}\hfill f_{X|Y}(x|y)dx& =\frac{f(x,y)dxdy}{f_Y(y)dy}\hfill \\ \hfill & \approx \frac{P\{x⩽X⩽x+dx,y⩽Y⩽y+dy\}}{P\{y⩽Y⩽y+dy\}}\hfill \\ \hfill & =P\{x⩽X⩽x+dx|y⩽Y⩽y+dy\}\hfill \end{array}$

In other words, for small values dx and dy, $f_{X|Y}(x|y)dx$ is approximately the conditional probability that X is between x and $x+dx$ given that Y is between y and $y+dy$ .

The conditional expectation of X, given that $Y=y$ , is defined for all values of y such that $f_Y(y)>0$ , by

$E[X|Y=y]={\int }_{-\infty }^{\infty }x{f}_{X|Y}(x|y)dx$

Example 3.5

Suppose the joint density of X and Y is given by

$f(x,y)=\{\begin{array}{cc}6xy(2-x-y),\hfill & 0<x<1,0<y<1\hfill \\ 0,\hfill & \text{otherwise}\hfill \end{array}$

Compute the conditional expectation of X given that $Y=y$ , where $0<y<1$ .

Solution:  We first compute the conditional density

$\begin{array}{cc}\hfill f_{X|Y}(x|y)& =\frac{f(x,y)}{f_Y(y)}\hfill \\ \hfill & =\frac{6xy(2-x-y)}{{\int }_0^{1}6xy(2-x-y)dx}\hfill \\ \hfill & =\frac{6xy(2-x-y)}{y(4-3y)}\hfill \\ \hfill & =\frac{6x(2-x-y)}{4-3y}\hfill \end{array}$

Hence,

$\begin{array}{cc}\hfill E[X|Y=y]& ={\int }_0^1\frac{6x^2(2-x-y)dx}{4-3y}\hfill \\ \hfill & =\frac{(2-y)2-\frac{6}{4}}{4-3y}\hfill \\ \hfill & =\frac{5-4y}{8-6y}■\hfill \end{array}$

Example 3.6 The t-Distribution

If Z and Y are independent, with Z having a standard normal distribution and Y having a chi-squared distribution with n degrees of freedom, then the random variable T defined by

$T=\frac{Z}{\sqrt{Y/n}}=\sqrt{n}\frac{Z}{\sqrt{Y}}$

is said to be a t-random variable with n degrees of freedom. To compute its density function, we first derive the conditional density of T given that $Y=y$ . Because Z and Y are independent, the conditional distribution of T given that $Y=y$ is the distribution of $\sqrt{n/y}Z$ , which is normal with mean 0 and variance $n/y$ . Hence, the conditional density function of T given that $Y=y$ is

$f_{T|Y}(t|y)=\frac{1}{\sqrt{2\pi n/y}}{e}^{-t^{2}y/2n}=\frac{y^{1/2}}{\sqrt{2\pi n}}{e}^{-t^{2}y/2n},-\infty <t<\infty$

Using the preceding, along with the following formula for the chi-squared density derived in Exercise 87 of Chapter 2,

$f_Y(y)=\frac{e^{-y/2}{y}^{n/2-1}}{2^{n/2}\mathrm{\Gamma }(n/2)},y>0$

we obtain the density function of T:

$f_T(t)={\int }_0^{\infty }{f}_{T,Y}(t,y)dy={\int }_0^{\infty }{f}_{T|Y}(t|y)f_Y(y)dy$

With

$K=\frac{1}{\sqrt{\pi n}{2}^{(n+1)/2}\mathrm{\Gamma }(n/2)},c=\frac{t^2+n}{2n}=\frac{1}{2}(1+\frac{t^2}{n})$

this yields

$\begin{array}{cc}\hfill f_T(t)& =\frac{1}{K}{\int }_0^{\infty }{e}^{-cy}{y}^{(n-1)/2}dy\hfill \\ \hfill & =\frac{c^{-(n+1)/2}}{K}{\int }_0^{\infty }{e}^{-x}{x}^{(n-1)/2}dx(\text{by letting}x=cy)\hfill \\ \hfill & =\frac{c^{-(n+1)/2}}{K}\mathrm{\Gamma }\left(\frac{n+1}{2}\right)\hfill \\ \hfill & =\frac{\mathrm{\Gamma }(\frac{n+1}{2})}{\sqrt{\pi n}\mathrm{\Gamma }(\frac{n}{2})}{(1+\frac{t^2}{n})}^{-(n+1)/2},-\infty <t<\infty ■\hfill \end{array}$

Example 3.7

The joint density of X and Y is given by

$f(x,y)=\{\begin{array}{cc}\frac{1}{2}y{e}^{-xy},\hfill & 0<x<\infty ,0<y<2\hfill \\ 0,\hfill & \text{otherwise}\hfill \end{array}$

What is $E[e^{X/2}|Y=1]$ ?

Solution:  The conditional density of X, given that $Y=1$ , is given by

$\begin{array}{cc}\hfill f_{X|Y}(x|1)& =\frac{f(x,1)}{f_Y(1)}\hfill \\ \hfill & =\frac{\frac{1}{2}{e}^{-x}}{{\int }_0^{\infty }\frac{1}{2}{e}^{-x}dx}=e^{-x}\hfill \end{array}$

Hence, by Proposition 2.1,

$\begin{array}{cc}\hfill E[e^{X/2}|Y=1]& ={\int }_0^{\infty }{e}^{x/2}{f}_{X|Y}(x|1)dx\hfill \\ \hfill & ={\int }_0^{\infty }{e}^{x/2}{e}^{-x}dx\hfill \\ \hfill & =2■\hfill \end{array}$

Example 3.8

Let $X_1$ and $X_2$ be independent exponential random variables with rates ${\mu }_1$ and ${\mu }_2$ . Find the conditional density of $X_1$ given that $X_1+X_2=t$ .

Solution:  To begin, let us first note that if $f(x,y)$ is the joint density of $X,Y$ , then the joint density of X and $X+Y$ is

$f_{X,X+Y}(x,t)=f(x,t-x)$

which is easily seen by noting that the Jacobian of the transformation

$g_1(x,y)=x,g_2(x,y)=x+y$

is equal to 1.

Applying the preceding to our example yields

$\begin{array}{cc}\hfill f_{X_1|X_1+X_2}(x|t)& =\frac{f_{X_1,X_1+X_2}(x,t)}{f_{X_1+X_2}(t)}\hfill \\ \hfill & =\frac{{\mu }_1{e}^{-{\mu }_{1}x}{\mu }_2{e}^{-{\mu }_2(t-x)}}{f_{X_1+X_2}(t)},0\le x\le t\hfill \\ \hfill & =C{e}^{-({\mu }_1-{\mu }_2)x},0\le x\le t\hfill \end{array}$

where

$C=\frac{{\mu }_1{\mu }_2{e}^{-{\mu }_{2}t}}{f_{X_1+X_2}(t)}$

Now, if ${\mu }_1={\mu }_2$ , then

$f_{X_1|X_1+X_2}(x|t)=C,0\le x\le t$

yielding that $C=1/t$ and that $X_1$ given $X_1+X_2=t$ is uniformly distributed on $(0,t)$ . On the other hand, if ${\mu }_1\ne {\mu }_2$ , then we use

$1={\int }_0^t{f}_{X_1|X_1+X_2}(x|t)dx=\frac{C}{{\mu }_1-{\mu }_2}(1-e^{-({\mu }_1-{\mu }_2)t})$

to obtain

$C=\frac{{\mu }_1-{\mu }_2}{1-e^{-({\mu }_1-{\mu }_2)t}}$

thus yielding the result:

$f_{X_1|X_1+X_2}(x|t)=\frac{({\mu }_1-{\mu }_2)e^{-({\mu }_1-{\mu }_2)x}}{1-e^{-({\mu }_1-{\mu }_2)t}}$

An interesting byproduct of our analysis is that

$f_{X_1+X_2}(t)=\frac{{\mu }_1{\mu }_2{e}^{-{\mu }_{2}t}}{C}=\{\begin{array}{cc}{\mu }^{2}t{e}^{-\mu t},\hfill & \text{if}{\mu }_1={\mu }_2=\mu \hfill \\ \frac{{\mu }_1{\mu }_2(e^{-{\mu }_{2}t}-e^{-{\mu }_{1}t})}{{\mu }_1-{\mu }_2},\hfill & \text{if}{\mu }_1\ne {\mu }_2\hfill \end{array}■$

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B9780128143469000081

Meta-Analysis and Latent Variable Models for Binary Data

Jian Qing Shi , in Handbook of Latent Variable and Related Models, 2007

2.3 Implementation

The basic idea is to treat the latent variable β = (β₁,…,β _k ) in (6) as missing and use an EM algorithm. Let Z be the collection of K vectors of z_i , the full log-likelihood for θ given ( Z, β) is

(10) $L\left(Z,\beta ;\theta \right)=\underset{i=1}{\overset{K}{\Sigma }}\left[log\left\{p\left(z_i|{\alpha }_i,{\beta }_i\right)\right\}+log\left\{\varphi \left({\beta }_i;{\mu }_{\beta },{\tau }_{\beta }\right)\right\}\right],$

where p(z_i | α _i , β _i ) is given by (5) but (α, β) is replaced by (α _i , β _i ). The EM algorithm involves the calculation of the expectation of the above full log-likelihood conditional on the current estimate of θ and the observation Z in the E-step, and then updates the estimate of θ by maximising this conditional expectation in the M-step. The details will be discussed in the rest of this section.

2.3.1 MCMC-EM: E step

In E-step of the (r + 1)th iteration, we need to calculate the following conditional expectation

(11) $\begin{array}{l}L\left(\theta |{\theta }^{\left(r\right)}\right)=E\left[L\left(Z,\beta ;\theta |Z,{\theta }^{\left(r\right)}\right)\right]\\ =\underset{i=1}{\overset{K}{\Sigma }}\left\{E\left[log\left(p\left(z_i|{\alpha }_i,{\beta }_i\right)\right)|Z,{\theta }^{\left(r\right)}\right]+E\left[log\left(\varphi \left({\beta }_i;{\mu }_{\beta },{\tau }_{\beta }\right)\right)|Z,{\theta }^{\left(r\right)}\right]\right\},\end{array}$

where θ ^(r) is the current estimate after the rth iteration. The expectation is calculated in terms of β. As there is no analytical form for Eq. (11), we use MCMC-EM algorithm; see, for example, Wei and Tanner (1990) and Booth and Hobert (1999).

To do this, we generate A number of random vectors { β ^a , a = 1,…,A} from the conditional distribution p(β| Z , θ ^(r)), and approximate (11) by

(12) $L_A\left(\theta |{\theta }^{\left(r\right)}\right)=\frac{1}{A}\underset{a=1}{\overset{A}{\Sigma }}L\left(Z,{\beta }^a;\theta \right).$

We will discuss how to update θ by maximising the above log-likelihood in Section 2.3.2. Now, we discuss how to generate those random numbers.

Given ( Z , θ), the latent variable β _i 's are conditional independent for i = 1,…, K. We can therefore generate random number for each component individually. For each component, its conditional density function is

(13) $p\left({\beta }_i|z_i,{\theta }^{\left(r\right)}\right)\propto p\left(z_i|{\alpha }_i,{\beta }_i\right)\varphi \left({\beta }_i;{\mu }_{\beta },{\tau }_{\beta }\right),$

where $p\left(z_i|{\alpha }_i^{\left(r\right)},{\beta }_i\right)$ is given by (5). We can use Metropolis-Hastings algorithm (Metropolis et al., 1953; Hastings, 1970) to generate a random variate β _i from the above density function. Carlin and Louis (2000) discussed the details for Metropolis-Hastings algorithm. Shi and Copas (2004) gave the details how to define a transition density and calculate acceptance probability for a similar problem.

2.3.2 MCMC-EM: M-step

In the M-step, we need to update θ by maximising the conditional likelihood (12). Since the unknown parameters (μ_β, τ_β) are involved in the second term only in the full log-likelihood function (10), the calculation of the maximum likelihood estimate is rather simple. This is to estimate (μ_β, τ_β) by maximising the following objective function

$\underset{i=1}{\overset{K}{\Sigma }}log\left\{\varphi \left({\beta }_i;{\mu }_{\beta },{\tau }_{\beta }\right)\right\}.$

It is actually equivalent to the problem that we have observed {β₁,…,β _K } from the normal distribution N(μ_β, τ_β) and want to calculate the maximum likelihood estimates of μ_β and τ_β. Their estimates are simply given by the sample mean and the sample standard error:

$\begin{array}{ccc}\overline{\beta }=\frac{1}{K}\Sigma {\beta }_i& \text{and}& V_{\beta }=\sqrt{\frac{1}{K}\Sigma {\left({\beta }_i-\overline{\beta }\right)}^2}.\end{array}$

Here, $\overline{\beta }$ and V_β are sufficient statistics for the parameters of (μ_β, τ_β). In the (r + 1)th iteration, M-step is to update β and τ_β by their conditional expectations $E\left(\overline{\beta }|Z,{\theta }^{\left(r\right)}\right)$ and $E\left(V_{\beta }|Z,{\theta }^{\left(r\right)}\right)$ . It is not possible to get an analytic form for those conditional expectations. We therefore use the random vectors {β ¹,…, β ^A } generated from the conditional distribution $p\left(\beta |Z,{\theta }^{\left(r\right)}\right)$ in E-step, and approximate these conditional expectations by their sample means:

$\begin{array}{ccc}\left({\overline{\beta }}^1+\cdots +{\overline{\beta }}^A\right)/A& \text{and}& \left(V_{\beta }^1+\cdots +V_{\beta }^A\right)\end{array}/A,$

where ${\overline{\beta }}^a$ and V_β ^a are the sample mean and sample standard error of ${\beta }^a=\left({\beta }_1^a,\dots ,{\beta }_K^a\right)$ , respectively.

The parameter α _i is updated by maximising $L_i=log\left\{p\left(z_i|{\alpha }_i,{\beta }_i\right)\right\}$ . There is no analytic solution. We use the following Newton method to approximate the estimate of α _i at the M-step. For simplifying the notation, we omit the index i here. We update α by the following subiteration

$\alpha ={\alpha }_0-\dot{L}\left({\alpha }_0\right)/\ddot{L}\left({\alpha }_0\right),$

where α₀ is the current estimate of α, $\dot{L}\left({\alpha }_0\right)$ and $\ddot{L}\left({\alpha }_0\right)$ are the first two derivatives of $L_i=log\left\{p\left(z_i|{\alpha }_i,{\beta }_i\right)\right\}$ in terms of α _i for the ith component. Bear in mind that we actually need to maximise the conditional expectation of $L=log\left\{p\left(z|\alpha ,\beta \right)\right\}$ given β. This can be approximated by the random numbers generated in MCMC-E step:

$L=\frac{1}{A}\underset{a=1}{\overset{A}{\Sigma }}log\left\{f\left(z|a,{\beta }^a\right)\right\}.$

The Newton method works very effectively as it is a univariate problem.

2.3.3 Average MCMC-EM algorithm and standard errors

As discussed in (Shi and Copas, 2002), the estimates by the MCMC-EM algorithm converges to the real maximum likelihood estimates when A is sufficiently large. The Monte Carlo error is mainly determined by the sample size A. To reduce the Monte Carlo error we need to take a sufficiently large A (see, e.g., Booth and Hobert, 1999). An alternative way is to use average MCMC-EM algorithm proposed in (Shi and Copas, 2002). Instead of increasing A, the average value of the estimates collected in the iterations after 'burn-in' is used. The average batch mean ${\overline{\theta }}^{\left(r\right)}$ is defined as the sample mean of $\left\{{\theta }^{\left(r-j+1\right)},{\theta }^{\left(r-j\right)},\dots ,{\theta }^{\left(r\right)}\right\}$ , where θ ^(r) is the estimate obtained in the rth iteration. The estimates ${\overline{\theta }}^{\left(r\right)}$ is roughly equivalent to the one calculated by using the MCMC-EM algorithm with Monte Carlo sample size J A, but the former is more efficient and much easier to implement than the latter.

The standard error of $\hat{\theta }$ can be calculated quite easily for the MCMC-EM algorithm. The related observed information matrix can be calculated by Louis (1982)

$I\left(\theta \right)=-E\left\{\ddot{L}\left(Z,\beta |Z\right)\right\}-E\left\{\dot{L}\left(Z,\beta |Z\right).{\dot{L}}^{\text{T}}\left(Z,\beta |Z\right)\right\},$

where $\dot{L}$ and $\ddot{L}$ are the first two derivatives of the full log-likelihood (10) with respect to θ. The expectation is defined in terms of β, which can be approximated by the random samples generated in MCMC-E step. In each iteration, we calculate

$I\left(\theta |\beta \right)=-\ddot{L}\left(Z,\beta \right)-\dot{L}\left(Z,\beta \right).{\dot{L}}^{\text{T}}\left(Z,\beta \right),$

evaluated at β = β ^a for each sample β ^a . The observed information matrix I (θ) is therefore approximated by the sample mean of $\left\{I\left(\theta |{\beta }^a\right),a=1,\dots ,A\right\}$ .

The final estimate of the information matrix can be calculated by the average 'batch mean' for average MCMC-EM algorithm. The variance of the estimates can be calculated from the inverse of the information matrix.

For the model we discussed in this section, we update the estimates of α _i 's and (μ_β, τ_β) in the M-step separately. Thus, it is rather simple to calculate standard errors for those parameters since they can also be calculated separately. The expressions for the first two derivatives related to α _i are (the index i is omitted for simplifying the notation)

$\begin{array}{l}\dot{L}\left(\alpha |\beta \right)=\underset{j=0}{\overset{n}{\Sigma }}\left[z_j-m_j\frac{exp\left(\alpha +\beta x_j\right)}{1+exp\left(\alpha +\beta x_j\right)}\right],\\ \ddot{L}\left(\alpha |\beta \right)=-\underset{j=0}{\overset{n}{\Sigma }}\left[m_j\frac{exp\left(\alpha +\beta x_j\right)}{{\left(1+exp\left(\alpha +\beta x_j\right)\right)}^2}\right],\end{array}$

with the related quantities evaluated for the ith component. The variance of α is therefore calculated by

${\left[\frac{1}{A}\underset{a=1}{\overset{A}{\Sigma }}\left\{-\ddot{L}\left(\alpha |{\beta }^a\right)-{\ddot{L}}^2\left(\alpha |{\beta }^a\right)\right\}\right]}^{-1}.$

The first two derivatives for (μ_β, τ_β) are given by

$\begin{array}{l}\dot{L}\left({\mu }_{\beta },{\tau }_{\beta }|\beta \right)=\left(\begin{array}{c}K\left(\overline{\beta }-{\mu }_{\beta }\right)/{\tau }_{\beta }^2\\ -K/{\tau }_{\beta }+{\Sigma }_{i=1}^K{\left({\beta }_i-{\mu }_{\beta }\right)}^2/{\tau }_{\beta }^3\end{array}\right),\\ \ddot{L}\left({\mu }_{\beta },{\tau }_{\beta }|\beta \right)=\left(\begin{array}{cc}-K/{\tau }_{\beta }^2& -2K\left(\overline{\beta }-{\mu }_{\beta }\right)/{\tau }_{\beta }^3\\ -2K\left(\overline{\beta }-{\mu }_{\beta }\right)/{\tau }_{\beta }^3& K/{\tau }_{\beta }^2-3{\Sigma }_{i=1}^m{\left({\beta }_i-{\mu }_{\beta }\right)}^2/{\tau }_{\beta }^4\end{array}\right),\end{array}$

where $\overline{\beta }$ is the average of {β₁,…,β _K }. The observed information matrix for (μ_β, τ_β) is calculated by

$\underset{a=1}{\overset{A}{\Sigma }}\left[-\ddot{L}\left({\mu }_{\beta },{\tau }_{\beta }|{\beta }^a\right)-\dot{L}\left({\mu }_{\beta },{\tau }_{\beta }|{\beta }^a\right){\dot{L}}^{\text{T}}\left({\mu }_{\beta },{\tau }_{\beta }|{\beta }^a\right)\right]/A.$

The final estimate of the information matrix can be approximated by the batch mean. The covariance matrix of (μ_β, τ_β) is the inverse of information matrix.

Read full chapter

URL:

https://www.sciencedirect.com/science/article/pii/B978044452044950015X

Goodness-of-fit test for interest rate models: An approach based on empirical processes

Abelardo Monsalve-Cobis , ... Manuel Febrero-Bande , in Computational Statistics & Data Analysis, 2011

2.1 Estimation of the model

There exists an extensive literature focused on estimating the drift and volatility functions of Model (2). For detailed reviews of the various methods, see Prakasa-Rao (1999) and Iacus (2008). In practise, all continuous time models are observed in discrete time. Therefore, to estimate of the parameters of the continuous time model (2), the discrete model (5) is considered. Given the Markovian nature of the process $r_t$ , which the discrete observations inherit, the likelihood function $L_n\left(\theta \right)$ of the discrete model can be derived from

(6) $L_n\left(\theta \right)=\prod _{i=0}^{n-1}{p}_{\theta }(\Delta ,r_{t_{i+1}}|r_{t_i})p_{\theta }\left(r_{t_0}\right)\text{,}$

where $p_{\theta }(\Delta ,r_{t_i+1}|r_{t_i})$ denotes the transition density or conditional density function. Usually $p_{\theta }$ , with $r_{t_0}=r_0$ , (initial condition). If the number of observations increases with time, it can be assumed that the relative weight of $p_{\theta }\left(r_{t_0}\right)$ in the likelihood function $L_n\left(\theta \right)$ decreases, and it can then be assumed that $p_{\theta }\left(r_{t_0}\right)=1$ . In the following we denote $r_0$ by $r_{t_0}$ . Let

(7) ${\ell }_n\left(\theta \right)=log{L}_n\left(\theta \right)=\sum _{i=0}^{n-1}log{p}_{\theta }(\Delta ,r_{t_{i+1}}|r_{t_i})+log\left(p_{\theta }\left(r_{t_0}\right)\right),=\sum _{i=0}^{n-1}{\ell }_i\left(\theta \right)+log\left(p_{\theta }\left(r_{t_0}\right)\right)\text{,}$

be the logarithm of the likelihood function, where ${\ell }_i\left(\theta \right)=log{p}_{\theta }(\Delta ,r_{t_{i+1}}|r_{t_i})$ , for $i=1,\dots ,n-1$ .

The maximum likelihood estimator (MLE) is expressed as $\hat{\theta }=arg{max}_{\theta }{\ell }_n\left(\theta \right)$ . If the parametric model that generates the observations ${\left\{r_{t_i}\right\}}_{{}_i}$ is known, then the maximum likelihood method is the natural choice. However, with the exception of particular cases, (e.g., the diffusion models of Vasicek (VAS)–Orsntein Uhlenbeck or CIR), an explicit form of the transition density function is not available. This is the case with the CKLS model, among others.

A technique commonly used to estimate Model (2) is to proceed as if the observations come from the Gaussian distribution, with a mean equal to the drift of the model and a standard deviation equal to the volatility function, and then to obtain the maximum likelihood estimator. Note that this method is efficient if $\Delta$ , the discretisation step, is sufficiently small (seeFan and Zhang, 2003), but the estimation can be significantly biased when $\Delta$ is large. To reduce the estimation bias, Aït-Sahalia (1999) suggests approximating the transition density function using Hermite polynomials. By applying this approximation, the author obtains a "root- $n$ " consistent maximum likelihood estimator of a diffusion process (seeAït-Sahalia, 2002).

This approximation is crucial to our goodness-of-fit method, and in the following, the maximum likelihood method suggested by Aït-Sahalia (2002) will be applied. The implementation of this estimator is quite reasonable and is in accordance with the requirements needed to develop the next section's proposal.

Read full article

URL:

https://www.sciencedirect.com/science/article/pii/S0167947311002052

goodethoughte.blogspot.com

Source: https://www.sciencedirect.com/topics/mathematics/conditional-density-function