The exact contents of citations.db main_text.text for this paper — one flattened LaTeX string, title through conclusion, appendix excluded, unmodified except for removing email addresses. This is what our citation measures are computed over.
98,398 characters
Exact Likelihood Inference and Robust Filtering for Gauss-Cauchy Convolution Models
\maketitle
\begin{abstract}
The convolution of a Gaussian and a Cauchy distribution, known as the Voigt distribution, is widely used in spectroscopy and provides a natural framework for modeling heavy-tailed measurement noise. We derive analytical expressions for its density, score, Hessian, and conditional moments using the scaled complementary error function, enabling stable maximum likelihood estimation without numerical convolution, finite-difference derivatives, or pseudo-Voigt approximations. The conditional expectation of the latent Gaussian component is governed by a redescending location score, so extreme observations are automatically discounted rather than propagated. This structure motivates the Gauss-Cauchy Convolution (GCC) filter for state-space models with Gaussian latent dynamics and heavy-tailed measurement errors. In an application to log realized volatility for the Technology Select Sector SPDR Fund, the GCC filter separates persistent latent variation from transient measurement noise and improves on Gaussian, Student-$t$, Huber, and related robust alternatives.
\end{abstract}
{\small\textit{{\noindent}Keywords:}}{\small{} Voigt profile, Convolutions, Heavy Tails, Robust Filtering, Kalman Filter.}{\small\par}
{\small\textit{{\noindent}JEL Classification:}}{\small{}
C13, C16, C22, C58, G17}{\small\par}
\setstretch{1.5}
\newpage
\section{Introduction}
The convolution of a Gaussian and a Cauchy distribution, known as the Voigt
distribution, arises naturally when finite-variance background noise is combined
with heavy-tailed disturbances. The Voigt profile is widely used in
spectroscopy, where Gaussian broadening is associated with Doppler effects and
Cauchy broadening with pressure or lifetime effects. It is also a natural
statistical model for measurement errors that combine ordinary noise with
occasional extreme observations.
Although an expression for the Voigt density has been known since \citet{Kendall_D:1938}, its use
as a likelihood-based statistical model has been limited by the perception that
exact inference is computationally difficult. Because the density is not an
elementary function, applications often rely on numerical convolution,
pseudo-Voigt approximations, nonlinear least squares, or simulation-based
methods. We show that these compromises are unnecessary for likelihood
inference. The scaled complementary error function provides a stable
representation of the Voigt density, and its differential identities yield
closed-form expressions for the likelihood derivatives.
This paper develops likelihood inference and filtering methods for the
Gauss-Cauchy convolution model. We first derive analytical expressions for the
density, score, Hessian, and conditional moments of the Voigt distribution.
These expressions reduce likelihood evaluation and differentiation to
evaluations of a standard special function and algebraic operations. We then
establish standard large-sample properties of the maximum likelihood estimator
and document its numerical performance in Monte Carlo simulations.
The same analytical structure is useful for robust signal extraction. If an
observation is the sum of a latent Gaussian signal and Cauchy measurement noise,
the conditional expectation of the Gaussian component is nonlinear in the
observation. The associated location score is approximately linear near the
center of the distribution but redescends in the tails. Hence, extreme
observations are automatically discounted rather than propagated into the latent
signal. This differs from the Kalman update, which is linear in the prediction error, and from Student-$t$ or Huber-type procedures, where redescending or bounded influence is obtained by specifying a robust observation density or loss
function directly. In the Gauss-Cauchy convolution, the redescending update is
instead an implication of signal extraction: it is the exact conditional mean of
the latent Gaussian component when the observation is contaminated by an additive Cauchy component. Thus, robustness is not imposed as an influence
function but derived from the convolution structure itself.
We use this result to construct the Gauss-Cauchy Convolution (GCC) filter for
linear state-space models with Gaussian latent dynamics and Voigt measurement
errors. The key step is that the Masreliez Gaussian prediction approximation
makes Tweedie's formula applicable to the state-prediction error, while the
Gaussian-Cauchy convolution structure keeps the prediction-error density inside
the Voigt family. Specifically, if the conditional state-prediction error is
approximated as Gaussian with variance $h_{t|t-1}$, then
$$
N(0,h_{t|t-1})+\mathcal V(0,\sigma,\gamma)
=
\mathcal V(0,\delta_t,\gamma),
\qquad
\delta_t^2=h_{t|t-1}+\sigma^2.
$$
Thus the state update is the Tweedie conditional-mean correction associated
with the Voigt prediction-error score. The filter nests the Kalman filter as
the special case $\gamma=0$ and the pure Cauchy measurement-error filter as the
special case $\sigma=0$. We assess the Masreliez approximation directly by
comparing the GCC recursion with an exact benchmark filter based on numerical
density propagation. The approximation error is small at both the
predictive-density level and the one-step-correction level, especially in the
empirically relevant region where the estimated Cauchy component is small
relative to the Gaussian component.
We illustrate the method using daily log realized volatility for the Technology Select Sector SPDR Fund (XLK). The realized-volatility series contains both persistent volatility movements and transient extreme observations associated with market stress, liquidity disruptions, and microstructure effects. The GCC filter separates these components by retaining a Gaussian core for ordinary measurement variation while using the Cauchy component to absorb occasional large deviations. In the empirical application, the GCC specification attains the largest pseudo-log likelihood among the competing filters considered, including Gaussian, Cauchy, Normal-Laplace, Student-$t$, and Huber alternatives.
The paper is related to several strands of literature. It contributes to work on convolution-based distributions and Voigt-profile estimation by showing that the special-function representation yields a tractable likelihood theory.
It is also related to robust filtering and score-driven dynamics
\citep{CrealKoopmanLucas:2013, Harvey:2013}, where heavy-tailed observation densities are used to reduce the influence of outliers. Recent contributions include robust score-driven filters based on Student-$t$ scores, including the multivariate filter of \citet{DInnocenzoLuatiMazzocchi:2023}, and approximate
filtering approaches connected to dynamic adaptive mixture models
\citep{CataniaDInnocenzoLuati:2026}. The GCC filter differs by deriving the redescending update as the conditional mean correction implied by the Gaussian-Cauchy convolution structure, rather than by specifying a heavy-tailed observation density or mixture structure directly.
The redescending conditional expectation also connects the GCC model to recent theoretical work on information aggregation with thin- and heavy-tailed signals. In global games, \citet{MorrisYildiz:2019} show that large shocks can trigger equilibrium shifts when agents rationally attribute extreme observations to a common heavy-tailed component rather than to idiosyncratic noise. A related ``too good to be true'' logic appears in \citet{WeibullMattssonVoorneveld:2007}, where sufficiently extreme signals are discounted because they are more likely to reflect noise than fundamentals. \citet{HautschHessMuller:2012} extend this exact mechanism to financial price discovery, demonstrating that fat-tailed noise causes rational agents to optimally discount extreme signals. The Gauss-Cauchy convolution provides a tractable parametric realization of this behavior: the conditional expectation of the latent Gaussian component is explicitly redescending, and the parameters $(\sigma,\gamma)$ dictate precisely when extreme observations are treated as valid signals and when they are discarded as noise.
The remainder of the paper is organized as follows. Section
\ref{sec:Voigt} defines the Gauss-Cauchy convolution and derives its conditional
moment structure. Section \ref{sec:MLE} develops likelihood inference for the
Voigt distribution. Section \ref{sec:GCCfilter} introduces the GCC filter.
Section \ref{sec:Simulations} evaluates the Masreliez approximation using exact
benchmark filtering. Section \ref{sec:Empirical} applies the filter to log
realized volatility. Section \ref{sec:Conclusion} concludes. Technical proofs
are collected in the Appendix.
\section{Gauss-Cauchy Convolution and Conditional Expectations}\label{sec:Voigt}
Let $Z\sim \mathcal{N}(0,\sigma^{2})$ and $X\sim\operatorname{Cauchy}(0,\gamma)$ be independent where $\sigma,\gamma>0$, and consider the convolution
$$
Y = \mu + Z + X.
$$
We denote the resulting distribution by $\mathcal{V}(\mu,\sigma,\gamma)$ because it is known as the \emph{Voigt profile} in the field
of spectroscopy, and we let $f_Y(y;\theta)$ denote its density, where $\theta=(\mu,\sigma,\gamma)^\prime$.
This density can be expressed in terms of the scaled complementary error function, $\mathsf{e}(w)
\equiv
\operatorname{erfcx}(w)$, where $\mathsf{e}(-iw)$ is the Faddeeva function. For $\operatorname{Re}(w)>0$, we have
$$
\mathsf{e}(w)
=
\frac{1}{\pi}
\int_{-\infty}^{\infty}
\frac{\exp(-t^2)}{w+it}dt,
$$
and for $Z\sim \mathcal{N}(0,1)$ we have
$\Pr(|Z|>r)=\mathsf{e}(\tfrac{r}{\sqrt{2}})e^{-r^2/2}$ for $r>0$.
We are particularly interested in evaluating $\mathsf{e}(w)$ on the vertical
\emph{Voigt line} in $\mathbb{C}$ traced out by $y\mapsto w_{y,\theta} = \tfrac{\gamma+i(y-\mu)}{\sigma\sqrt{2}}$, and we introduce the notation for the real and imaginary parts of $\mathsf{e}(w_{y,\theta})$
$$
\mathsf{u}(y;\theta)
\equiv
\operatorname{Re}[\mathsf{e}(w_{y,\theta})],
\qquad
\mathsf{v}(y;\theta)
\equiv
\operatorname{Im}[\mathsf{e}(w_{y,\theta})],
$$
so that $\mathsf{e}(w_{y,\theta})=
\mathsf{u}(y;\theta)+i\mathsf{v}(y;\theta)$.
\begin{proposition}[Kendall, 1938]
\label{prop:Voigt}The density of the Gauss-Cauchy convolution, $\mathcal{V}(\mu,\sigma,\gamma)$, can be expressed
\begin{equation}
f_Y(y;\theta)=\frac{1}{\sqrt{2\pi\sigma^{2}}}u(y;\theta),\qquad \theta=(\mu,\sigma,\gamma)^\prime.\label{eq:Kendall}
\end{equation}
\end{proposition}
The expression for this density was first obtained by \citet{Kendall_D:1938}. We include a proof of Proposition \ref{prop:Voigt} in the appendix based on the characteristic function $\varphi_{Z+X}(s)=e^{-s^{2}/2-|s|}$. Figure \ref{fig:DenCauchyNormal} shows the density of $Y$ for $\mu=0$ and $\sigma=\gamma=1$, together with the best approximating Student-$t$ distribution, as defined by the Kullback-Leibler discrepancy. The visible differences highlight the specification error from using Student-$t$ distributions as proxies.\footnote{The Voigt density can be viewed as a limiting case of the convolution of two Student-$t$ distributions, a class of convolution problems known to be analytically difficult. For recent progress in this area, see \citet{Nason:2006}, \citet{Forchini:2008}, \citet{BergVignat:2010}, and \citet{HansenTong:2026}.} A second common approximation is the pseudo-Voigt profile, a two-component mixture of a Gaussian density and a Cauchy density with the same location.\footnote{The pseudo-Voigt density takes the form $(1-\eta_p)\phi(y;\mu,\sigma_p)+\eta_p c(y;\mu,\gamma_p)$, with $0\leq\eta_p\leq1$, where $\phi$ is the Gaussian density and $c$ is the Cauchy density.} It provides a closer density approximation in this example, but it is not a convolution and therefore does not preserve the additive-noise interpretation used below.
The Voigt profile uniquely captures the superposition of finite-variance noise and infinite-variance jumps, a geometry that generic heavy-tailed approximations fail to reproduce.
\begin{figure}[tbh]
\begin{centering}
\includegraphics[width=0.49\textwidth]{Figures/VoigtWithApproxs.pdf}
\includegraphics[width=0.5\textwidth]{Figures/VoigtApproxTails.pdf}
\par\end{centering}
\caption{{\small The left panel shows the Voigt density, $\mathcal{V}(0,1,1)$, together with the best approximating Student-$t$ distribution, $(\sigma_t,\nu_t)\approx(1.45,1.22)$, and pseudo-Voigt mixture, $(\eta_p,\sigma_p,\gamma_p)\approx(0.65,1.62,1.65)$. The right panel shows the corresponding Q-Q plots. The approximation parameters are chosen to minimize the Kullback-Leibler discrepancy, but visible discrepancies remain.\label{fig:DenCauchyNormal}}}
\end{figure}
Interestingly, the density (\ref{eq:Kendall}) can be expressed from the Mills ratio for a standard normal random variable, $m(t)=\Phi(-t)/\phi(t)$,
where $\phi$ and $\Phi$ are the density and cumulative distribution function for $\mathcal{N}(0,1)$, respectively.\footnote{We thank George Tauchen for pointing out the connection between the Voigt representation and the Mills ratio.} This representation highlights
the connection between the Voigt profile and the hazard rate of the normal distribution, providing another avenue for numerical evaluation.
\begin{corollary}\label{cor:MillsRatio}
The Voigt density can be expressed in terms of the complex Mills ratio for a standard normal variable:
$$
f_{Y}(y;\theta)=\frac{1}{\pi\sigma}\operatorname{Re}\left[m(\tfrac{\gamma+i(y-\mu)}{\sigma})\right],\qquad \theta=(\mu,\sigma,\gamma)^\prime,
$$
where $f_{Y}$ is the density in (\ref{eq:Kendall}), and $m(t)=\Phi(-t)/\phi(t)$.
\end{corollary}
\begin{lemma}\label{lem:erfcx}
Let $\sigma>0$ and $\gamma>0$. Then: (i) $\mathsf{e}$ is entire and satisfies
$$
\mathsf{e}^\prime(w) = 2w\mathsf{e}(w)-\frac{2}{\sqrt{\pi}}.
$$
Moreover, for every integer $n\geq0$, there exist polynomials $p_n$ and
$q_n$ such that
$$
\mathsf{e}^{(n)}(w)=p_n(w)\mathsf{e}(w)+q_n(w).
$$
(ii) The real and imaginary parts satisfy the symmetry relations
$\mathsf{u}(\mu+r;\theta)=\mathsf{u}(\mu-r;\theta)$,
$\mathsf{v}(\mu+r;\theta)=-\mathsf{v}(\mu-r;\theta)$,
and $\mathsf{u}(y;\theta)>0$ for all real $y$.
\end{lemma}
Lemma \ref{lem:erfcx} collects the special-function identities that make the Voigt likelihood tractable. Part (i) states that $\mathsf{e}=\operatorname{erfcx}$ is closed under differentiation in the sense that every derivative is an algebraic function of $w$ and $\mathsf{e}(w)$. Part (ii) gives the positivity and symmetry properties needed for the real and imaginary components along the Voigt line.
The positivity of $\mathsf{u}(y;\theta)$ ensures that the density, log-density, score, Hessian, and conditional moment formulas involving ratios such as
$$
\frac{\mathsf{v}(y;\theta)}{\mathsf{u}(y;\theta)}
$$
are well defined for all real $y$. The symmetry relations explain the corresponding even-odd structure of the density and the conditional moments derived below.
The main likelihood implication of Lemma \ref{lem:erfcx} is the following algebraic closure property.
\begin{corollary}[Derivatives of Voigt density]\label{cor:VoigtAlgebraicClosure}
All finite-order derivatives of the Voigt density,
$f_Y(y;\theta)$, with respect to $y$, $\mu$, $\sigma$, and $\gamma$ are
algebraic functions of $w_{y,\theta}$, $\mathsf{u}(y;\theta)$,
$\mathsf{v}(y;\theta)$, and the parameters $\theta=(\mu,\sigma,\gamma)^\prime$. Consequently, all finite-order derivatives of the log-density $\log f_Y(y;\theta)$ can be expressed as rational functions of the same quantities, with denominators involving only powers of $\mathsf{u}(y;\theta)$ and $\sigma$. Since $\mathsf{u}(y;\theta)>0$ and $\sigma>0$, these expressions are well defined for all real $y$.
\end{corollary}
Corollary \ref{cor:VoigtAlgebraicClosure} is the central tractability result for likelihood inference. Although the Voigt density is defined as a Gaussian-Cauchy convolution, the scaled complementary error function absorbs the otherwise unstable exponential term and turns differentiation of the density into algebraic operations on a single complex-valued function. Thus the likelihood, score, Hessian, and higher-order derivatives can all be evaluated from the same pair of real-valued quantities, $\mathsf{u}(y;\theta)$ and $\mathsf{v}(y;\theta)$, without numerical convolution, quadrature, finite-difference approximation, or pseudo-Voigt approximations. This is why exact maximum likelihood estimation is more direct than the convolution representation might suggest.
\subsection{Conditional Distribution, Expectation, and Variance}
We now consider the signal-extraction problem implied by the convolution
$$
Y=\mu+Z+X.
$$
The objective is to infer the latent Gaussian component $Z$ from the observed
realization of $Y$, treating the Cauchy component $X$ as measurement noise.
This formulation is directly relevant for robust filtering, where $Z$
represents the latent signal of interest and $X$ captures transient
heavy-tailed contamination. The following theorem gives the closed form conditional
density of $Z$ given $Y=y$.
\begin{theorem}
\label{thm:CondDensity}Suppose that $Z\sim \mathcal{N}(0,\sigma^{2})$ and
$X\sim\text{\ensuremath{\mathrm{Cauchy}(0,\gamma)}}$ are independent
and define $Y=\mu+Z+X$. Then the conditional density of $Z$ given $Y=y$
is given by
$$
f_{Z|Y}(z|y)=\frac{1}{\gamma\pi }\frac{1}{\mathsf{u}(y;\theta)}\frac{\exp(-\frac{1}{2}z^{2}/\sigma^{2})}{1+(\frac{y-\mu-z}{\gamma})^{2}},
$$
where $f_{Z|Y}(z|y)\rightarrow f_{Z}(z)=e^{-z^{2}/(2\sigma^{2})}/\sqrt{2\pi\sigma^{2}}$
as $(y-\mu)\rightarrow\pm\infty$.
\end{theorem}
This conditional density is shown in Figure \ref{fig:Conditional-density}
for selected values of $y$ with $\sigma=\gamma=1$. The vertical
lines along the $x$-axis indicate the conditional expected value
for each conditional density.
For the special case, $y=\mu$, the conditional density is
$$
f_{Z|Y}(z|y=\mu)=\frac{1}{\gamma\pi}\frac{1}{\Pr(|Z|\geq\tfrac{\gamma}{\sigma})}\frac{\exp(-\frac{1}{2}\frac{z^{2}+\gamma^{2}}{\sigma^{2}})}{1+(\frac{z}{\gamma})^{2}}.
$$
\begin{figure}[t]
\begin{centering}
\includegraphics[width=0.8\textwidth]{Figures/CondDensities}
\par\end{centering}
\caption{{\small Conditional density of $Z$ given $Y=y$, for selected values
of $y$, where $Y=Z+X$ with $Z\sim \mathcal{N}(0,1)$ and $X\sim\mathrm{Cauchy(0,1)}$
independent. The corresponding conditional expectations, $\mathbb{E}[Z|Y=y]$,
are indicated with vertical lines. The unconditional density of $Z$
(standard normal) is represented with the shaded area. }\label{fig:Conditional-density}}
\end{figure}
The conditional moments of $Z$ given $Y$ can be obtained from a Gaussian-convolution version of Tweedie's formula. The key point is that the conditional mean of an additive Gaussian component is determined by the score of the marginal density, and the conditional variance by its derivative. See \citet{Efron:2011} for a modern treatment of Tweedie's formula and \citet{Robbins:1956} for the original published attribution.
\begin{proposition}[Tweedie's formula for Gaussian convolutions]\label{prop:TweedieGaussian}
Let $Z\sim\mathcal{N}(0,\sigma^2)$, with $\sigma>0$, be independent of a proper real-valued random variable $X$, and define $Y=\mu+Z+X$.
Then $Y$ has a strictly positive and smooth density $f_Y$. Moreover, for all real $y$,
$$
\mathbb{E}[Z|Y=y] = -\sigma^2\frac{\partial}{\partial y}\log f_Y(y),
$$
and
$$
\operatorname{var}(Z|Y=y) = \sigma^2+\sigma^4\frac{\partial^2}{\partial y^2}\log f_Y(y) = \sigma^2\left(1-\frac{\partial}{\partial y}\mathbb{E}[Z|Y=y]\right).
$$
\end{proposition}
In the Gauss-Cauchy convolution, Proposition~\ref{prop:Voigt} gives the marginal density $f_Y$ explicitly in terms of $\mathsf{u}(y;\theta)$. Therefore Proposition~\ref{prop:TweedieGaussian} reduces the conditional-moment calculation to differentiating the Voigt log-density. Using the score and Hessian identities in Lemma~\ref{lem:ScoreHessian} gives the following closed-form expressions.
\begin{corollary}
\label{cor:CondExpect}In the Gauss-Cauchy convolution model, the conditional expectation and conditional variance of the latent Gaussian component $Z$ given $Y=y$ are
\begin{eqnarray}
\mathbb{E}\left[Z|Y=y\right]&=&(y-\mu) +\gamma\frac{\mathsf{v}(y;\theta)}{\mathsf{u}(y;\theta)}\label{eq:CondExpect}\\
\mathbb{V}(Z|Y=y) &=& \sqrt{2/\pi}\frac{\sigma\gamma}{\mathsf{u}(y;\theta)} - \gamma^2\left(1 + \frac{\mathsf{v}^{2}(y;\theta)}{\mathsf{u}^{2}(y;\theta)}\right)
\end{eqnarray}
and the conditional expectation reaches its maximum and minimum values at $y^{*}$ satisfying
$$
(\gamma^{2}+\sigma^{2})\mathsf{u}^{2}(y^{\ast};\theta)+\gamma^{2}\mathsf{v}^{2}(y^{\ast};\theta)-\sqrt{2/\pi}\sigma\gamma \mathsf{u}(y^{\ast};\theta)=0.
$$
\end{corollary}
The conditional expectation of $Z$ given
$Y$ is antisymmetric about $\mu$ and, conversely, the conditional variance is symmetric about $\mu$, driven by the fact that it depends on the prediction error strictly through the even function $\mathsf{u}(y)$ and the squared odd function $\mathsf{v}(y)^2$.
The conditional expectation and variance are plotted for $(\mu,\gamma,\sigma)=(0,1,1)$
in Figure \ref{fig:Conditional-expectation}. For this case, the maximum conditional expectation is $\max_{y}\mathbb{E}\left[Z|Y=y\right]\approx0.7486$
for $y\approx2.4637$ and the conditional variance ranges between $\mathbb{V}(Z|Y=0)\approx 0.5251$ and $\mathbb{V}(Z|Y=y)\approx 1.1603$ for $y^\ast \approx \pm 3.6621$.\footnote{That
$\operatorname{var}(Z|Y=y)$ can exceed $\operatorname{var}(Z)$ is not inconsistent with the law of total variance, because the latter implies only that
$\mathbb{E}[\operatorname{var}(Z|Y)]\leq\operatorname{var}(Z)$, not that
$\operatorname{var}(Z|Y=y)\leq\operatorname{var}(Z)$ for every $y$. In the standardized case, $\sigma=\gamma=1$ and $\mu=0$, the thresholds for maximum signal extraction ($y\approx2.4637$) and maximum signal confusion ($y\approx3.6621$) are the positive roots of the second and third derivatives of the Voigt log-likelihood, respectively.}
As $y\rightarrow\pm\infty$
the conditional moments converge to the unconditional moments,
$\mathbb{E}\left[Z|Y=y\right]\rightarrow 0$
and $\mathbb{V}\left[Z|Y=y\right]\rightarrow 1$. The reason is that the conditional distribution
converges to the unconditional distribution as $y$ increases in absolute value, i.e. $f_{Z|Y}(z|y)\rightarrow f_{Z}(z)$ (almost everywhere) as $y\rightarrow\pm\infty$.
The intuition for this is simply that the tail of a Gaussian vanishes much faster than that of a Cauchy, such that a very extreme realization of $y$ is almost fully attributed to the Cauchy component.
Figure \ref{fig:Conditional-expectation} provides deep insight into the filter we propose later. The shape of (\ref{eq:CondExpect}) is approximately linear near the origin similar to a standard Kalman filter. However, for large values of $|y|$ the filter enters a nonlinear regime, and the expectation does not merely saturate but redescends toward zero. This non-monotonicity is key and implies that the filter identifies observations that are ``too large to be true signals'' and isolates them as pure noise, protecting the latent trend estimate from contamination.
Since $Y=\mu + Z+X$ we may define
$\mathbb{E}[X|Y]= -\gamma \frac{\mathsf{v}(Y)}{\mathsf{u}(Y)}$, despite $X$ being Cauchy distributed, with $\mathbb{E}|X|=\infty$.
The conditional mean in Figure \ref{fig:Conditional-expectation} bears a striking resemblance to the score function of the Student-$t$ distribution shown in \citet[Figure 1]{HarveyLuati:2014}. This is noteworthy because the two objects arise from different constructions. In the score-driven literature, the redescending shape is induced by specifying a heavy-tailed observation density, such as the Student-$t$, whose score limits the influence of large observations. In contrast, the redescending shape in Figure \ref{fig:Conditional-expectation} emerges endogenously as the exact conditional mean of a latent Gaussian signal observed with Cauchy noise. Thus, the Gauss-Cauchy convolution provides a probabilistic micro-foundation for robust update functions of the type used in score-driven models.
\begin{figure}[tbh]
\begin{centering}
\includegraphics[width=0.95\textwidth]{Figures/CondMoments.pdf}
\par\end{centering}
\caption{{\small The left panel shows the conditional expectation of $Z$ given $Y=y$, where $Y=Z+X$ with $Z\sim\mathcal N(0,1)$ and $X\sim\operatorname{Cauchy}(0,1)$ independent. The right panel shows the corresponding conditional variance. The insets present the same quantities over a wider range of values of $y$.\label{fig:Conditional-expectation}}}
\end{figure}
The non-monotonicity of the conditional variance is inextricably linked to the redescending nature of the conditional expectation. By differentiating Tweedie's formula for the conditional mean, $\mathbb{E}[Z | Y=y] = -\sigma^2 \frac{\partial \log f_Y(y)}{\partial y}$, with respect to $y$, we obtain $\frac{\partial}{\partial y}\mathbb{E}[Z| Y=y] = -\sigma^2 \frac{\partial^2 \log f_Y(y)}{\partial y^2}$. Substituting this into the conditional variance equation yields the exact identity:
\begin{equation}
\mathbb{V}(Z| Y=y) = \sigma^2 \left( 1 - \frac{\partial}{\partial y} \mathbb{E}[Z| Y=y] \right).
\end{equation}
This establishes three distinct filtering regimes based on the slope of the conditional expectation. Near zero, where $\mathbb{E}(Z|Y=y)$ is increasing in $|Y|$, the observation is deemed reliable, and the conditional variance is strictly less than the unconditional variance $\sigma^2$, reflecting a reduction in uncertainty.
Then, at the exact local extrema of the conditional expectation ($y \approx \pm 2.4637$ for the standard case), the slope is zero, and the conditional variance exactly equals $\sigma^2$. This represents the exact threshold where the filter ceases to trust the magnitude of the innovation. As the observation moves further into the tails, the expectation redescends and the conditional variance strictly exceeds $\sigma^2$. The increased uncertainty occurs because the filter cannot confidently attribute the moderate outlier to either the Gaussian signal or the Cauchy noise. As $|y| \to \infty$, $\frac{\partial}{\partial y}\mathbb{E}(Z|Y=y) \to 0$ and the conditional variance smoothly asymptotes back to the unconditional variance as the observation is entirely discounted.
\subsection{Higher-Order Conditional Cumulants}
Higher-order conditional cumulants can be obtained in the same way. The first step is not specific to the Voigt model, but holds for any Gaussian convolution. Let
$$
K(s|y)=\log\mathbb{E}\left[\exp(sZ)|Y=y\right]
$$
denote the conditional cumulant generating function of $Z|Y=y$. If $Y=\mu+Z+X$, where $Z\sim\mathcal{N}(0,\sigma^2)$ is independent of an arbitrary proper random variable $X$, then completing the square gives
$$
K(s|y)=\frac{\sigma^2s^2}{2}+\log f_Y(y-\sigma^2s)-\log f_Y(y).
$$
The conditional cumulants are defined by
$$
\kappa_r(y)=\left.\frac{\partial^r}{\partial s^r}K(s|y)\right|_{s=0}.
$$
Therefore
$$
\kappa_1(y)=-\sigma^2\frac{\partial}{\partial y}\log f_Y(y),
\qquad
\kappa_2(y)=\sigma^2+\sigma^4\frac{\partial^2}{\partial y^2}\log f_Y(y).
$$
These are the conditional mean and variance identities in Proposition \ref{prop:TweedieGaussian}. For $r\geq3$,
$$
\kappa_r(y)=(-\sigma^2)^r\frac{\partial^r}{\partial y^r}\log f_Y(y),
$$
and equivalently,
$
\kappa_{r+1}(y)=-\sigma^2\frac{\partial}{\partial y}\kappa_r(y)$, for all $r\geq2$.
For the Voigt Gauss-Cauchy convolution, Proposition \ref{prop:Voigt} gives $f_Y$ explicitly in terms of $\mathsf{u}(y;\theta)$. Hence all higher-order conditional cumulants are obtained by differentiating the closed-form variance expression in Corollary \ref{cor:CondExpect}. By Lemma \ref{lem:erfcx}, these derivatives are again algebraic functions of $\mathsf{u}(y;\theta)$, $\mathsf{v}(y;\theta)$, and the parameters.
The tractability is therefore not limited to the conditional mean and variance. Although the Voigt random variable $Y$ has Cauchy tails and hence does not have finite positive integer moments when $\gamma>0$, the conditional distribution of the latent Gaussian component $Z|Y=y$ has moments of all orders. Thus conditional skewness, kurtosis, and higher conditional moments of the latent Gaussian component can be evaluated without numerical integration.
\section{Likelihood-Based Inference for Voigt Distributions}\label{sec:MLE}
Having established the exact density and conditional properties of the Voigt distribution, we now turn to the problem of parameter estimation. While the Voigt distribution is fundamental in spectroscopy and physics, rigorous statistical inference for its parameters has historically been underdeveloped. The closed-form expression in (\ref{eq:Kendall}) allows us to proceed directly with Maximum Likelihood Estimation (MLE), avoiding the approximations commonly used in the existing literature. These results are obtained by leveraging the differential properties of the scaled complementary error function, which leads to exact analytical expressions for the score vector and Hessian matrix.
\begin{lemma}
\label{lem:ScoreHessian}Let $Y\sim\operatorname{Voigt}(\mu,\sigma,\gamma)$,
then the score and Hessian matrix are given by
\begin{align*}
s_{\mu} & =\tfrac{1}{\sigma^2}\left(y-\mu+\gamma\tfrac{v}{u}\right)=\tfrac{1}{\sigma^2}\mathbb{E}(Z|Y=y)\\
s_{\sigma} & =\tfrac{1}{\sigma^{3}u}\left((\tilde{y}^{2}-\gamma^{2}-\sigma^{2})u+2\gamma\tilde{y}v+\sqrt{\tfrac{2}{\pi}}\sigma\gamma\right)\\
s_{\gamma} & =\tfrac{1}{\sigma^{2}u}\left(\gamma u-\tilde{y}v-\sqrt{\tfrac{2}{\pi}}\sigma\right),
\end{align*}
where $\tilde{y}=y-\mu$, $u=\mathsf{u}(y;\mu,\sigma,\gamma)$ and $v=\mathsf{v}(y;\mu,\sigma,\gamma)$ and the elements of the Hessian,
$ H_\theta $, are
\begin{align*}
H_{\mu\mu} &= \tfrac{1}{\sigma}s_{\sigma}-s_{\mu}^{2}&\quad&
H_{\mu\sigma} =-\tfrac{1}{\sigma}\left(s_{\mu}+\gamma H_{\mu\gamma}-\tilde{y}H_{\mu\mu}\right) \\
H_{\gamma\gamma}& = {-\tfrac{1}{\sigma}s_{\sigma}-s_{\gamma}^{2}}&\quad&
H_{\gamma\sigma} = -\tfrac{1}{\sigma}\left(s_{\gamma}+\gamma H_{\gamma\gamma}-\tilde{y}H_{\mu\gamma}\right)\\
H_{\mu\gamma} &= \tfrac{1}{\sigma^{2}}\left(\tilde{y}s_{\gamma}+\gamma s_{\mu}+\tfrac{v}{u}\right)-s_{\mu}s_{\gamma}&\quad&
H_{\sigma\sigma} = -\tfrac{1}{\sigma}\left(s_{\sigma}+\gamma H_{\gamma\sigma}-\tilde{y}H_{\mu\sigma}\right).
\end{align*}
\end{lemma}
One significance of Lemma \ref{lem:ScoreHessian} is computational. The
perceived intractability of the Voigt integral has led researchers to rely on
computationally expensive methods or approximation-based substitutes.\footnote{For example, \citet{CannasPiras:2025} note that direct numerical integration is time consuming and discuss the use of pseudo-Voigt approximations and Bayesian methods. The Voigt profile also has a substantial analytic literature on special-function representations, including Kummer-function expansions \citep{DiRoccoAguirreTellez:2004} and Mellin-Barnes, Fox $H$-, and Meijer $G$-function representations \citep{PagniniSaxena:2008}. Related computational bottlenecks have also motivated workarounds in physics and spectroscopy, including deep learning approaches \citep{StemockChurchillLeeHassanDoughtyOchoa:2024,
AliVanZijlPrasuhnWirestamKnutssonYadav:2025}. Lemma \ref{lem:ScoreHessian}
shows that these complications are not intrinsic to likelihood-based inference
for the Voigt model.} However, once the density is represented through the scaled
complementary error function, the likelihood, score, and Hessian are all
computed from the same quantities $\mathsf{u}(y;\theta)$ and
$\mathsf{v}(y;\theta)$. Exact maximum likelihood estimation therefore reduces
to repeated evaluation of a standard special function and algebraic operations,
without numerical convolution, finite-difference derivatives, or pseudo-Voigt
approximations.
A natural concern is whether standard likelihood-based inference remains valid
for a model whose distribution has Cauchy tails. Since the Cauchy component dominates the tails of $Y$, the distribution has no finite positive integer moments; in particular, even the mean is undefined. One might therefore expect conventional asymptotic normality to fail. The relevant issue for maximum likelihood estimation, however, is not whether $Y$ has finite moments, but whether the score, Hessian, and Fisher information are well defined. In the Gauss-Cauchy convolution model they are, so exact maximum likelihood remains a regular estimation problem despite the absence of moments for $Y$. This is analogous to cases that arise in GARCH models, where asymptotic normality can hold even when the observed process has limited or no finite moments; see \citet{JensenRahbek:2004}. Theorem \ref{thm:MLE-consistent-asN} establishes consistency and asymptotic normality of the MLE for the Voigt model.
\begin{assumption}[Compactness]\label{assu:Compact}
The parameter vector satisfies $\theta\in\Theta$, where $\Theta$ is the compact set
$$
\Theta
=
\left\{
(\mu,\sigma,\gamma)\in\mathbb{R}^{3}:
|\mu|\leq\mu_{\max},
0<\sigma_{\min}\leq\sigma\leq\sigma_{\max},
0<\gamma_{\min}\leq\gamma\leq\gamma_{\max}
\right\}.
$$
\end{assumption}
Assumption \ref{assu:Compact} restricts attention to the regular Voigt case in which both scale parameters are strictly positive and bounded away from zero. This excludes the boundary cases $\sigma=0$ and $\gamma=0$, corresponding to the pure Cauchy and pure Gaussian limits, respectively. These limiting cases are mathematically meaningful, but they require separate treatment because the analytical expressions below involve ratios and derivatives evaluated under $\sigma>0$ and $\gamma>0$.
The lower bound on $\gamma$ is also useful for the uniform integrability condition used in the consistency proof. If the true distribution has Cauchy tails but the parameter space includes the pure Gaussian boundary $\gamma=0$, the Gaussian tail can make $|\log f_Y(Y;\theta)|$ grow quadratically in $Y$ for some admissible parameter values, while $Y$ itself has Cauchy tails. Bounding $\gamma$ away from zero avoids this difficulty and yields a common logarithmic envelope for the likelihood. The lower bound on $\sigma$ similarly keeps the analysis within the interior of the regular Voigt family, where the score and Hessian formulas are ordinary derivatives with respect to $(\mu,\sigma,\gamma)$.
\begin{theorem}[Consistency and Asymptotic Normality of the Voigt MLE]
\label{thm:MLE-consistent-asN}
Let $Y_1,\ldots,Y_n$ be iid with distribution
$\mathcal{V}(\mu_0,\sigma_0,\gamma_0)$, and let
$\theta_0=(\mu_0,\sigma_0,\gamma_0)^\prime\in\Theta$, where $\Theta$ is the
compact parameter space defined in Assumption \ref{assu:Compact}. Define the
maximum likelihood estimator by
$$
\hat\theta_n
=
\arg\max_{\theta\in\Theta}\ell_n(\theta),
\qquad
\ell_n(\theta)
=
\frac{1}{n}\sum_{i=1}^n\log f_Y(Y_i;\theta).
$$
Then $\hat\theta_n\xrightarrow{p}\theta_0$.
If, in addition, $\theta_0$ is an interior point of $\Theta$, then the
information matrix equality holds:
$$
\mathcal{I}_{\theta_0} = \mathbb{E}_{\theta_0}\left[s(Y;\theta_0)s(Y;\theta_0)^\prime\right]
= -\mathbb{E}_{\theta_0}\left[H(Y;\theta_0)\right],
$$
where
$$ s(y;\theta)=
\frac{\partial\log f_Y(y;\theta)}{\partial\theta},
\qquad H(y;\theta)
= \frac{\partial^2\log f_Y(y;\theta)}
{\partial\theta\partial\theta^\prime}.
$$
Moreover, $\mathcal{I}_{\theta_0}$ is nonsingular and
$\sqrt{n}(\hat\theta_n-\theta_0)
\xrightarrow{d}
\mathcal{N}(0,\mathcal{I}_{\theta_0}^{-1})$.
Finally, by symmetry of the Voigt density, $\mathcal{I}_{\theta_0}$ is
block diagonal with
$\mathcal{I}_{\mu_0\sigma_0}=\mathcal{I}_{\mu_0\gamma_0}=0$.
\end{theorem}
Despite the lack of moments, $\mathbb{E}|Y|=\infty$,
Theorem \ref{thm:MLE-consistent-asN} shows that likelihood-based inference is quite standard in the Gauss-Cauchy convolution model. Convergence in distribution is at the conventional root-$n$ rate and the limit distribution is Gaussian, such that standard errors and $t$-statistics have their usual interpretation. The underlying reason is that the likelihood function and its derivatives are well-behaved. Moreover, the information-matrix equality holds, and the simplest way to estimate $\mathcal{I}_{\theta_{0}}$ appears to be to evaluate $\mathbb{E}_{\theta_{0}}[s_{\theta_{0}}s_{\theta_{0}}^{\prime}]$ by numerical integration.
The asymptotic independence of $\hat{\mu}$ and $(\hat{\sigma},\hat{\gamma})$ is a consequence of $s_\mu$ being an odd function of the $\tilde{y}=y-\mu$, whereas $s_\gamma$ and $s_\sigma$ are both even functions. Consequently, ($\sigma, \gamma$) will have the same asymptotic variances in the simpler model where $\mu$ is known.
Moreover, the asymptotic covariance matrix for $(\sigma,\gamma)$ is
homogeneous of degree two. Hence the ratio of their asymptotic standard
deviations depends only on $\lambda=\gamma/\sigma$. We define this ratio
as\footnote{The ratio $\lambda$ is proportional to the Voigt damping parameter
used in astrophysics, $a=\lambda/\sqrt{2}$.}
$$
R(\lambda)\equiv
\frac{\operatorname{aStd}(\sigma;\mu,\sigma,\lambda\sigma)}
{\operatorname{aStd}(\gamma;\mu,\sigma,\lambda\sigma)}.
$$
Figure \ref{fig:seratio} reports the resulting asymptotic standard deviations
and their ratio. The left panel shows the asymptotic standard deviations for
$\sigma$ and $\gamma$ as functions of $\lambda=\gamma/\sigma$, under the
normalization $(\sigma,\gamma)=(1,\lambda)$. Both axes are log transformed. The
relationship is close to log-linear for $\gamma$, corresponding to an
approximately constant elasticity, whereas the relationship for $\sigma$ is
convex.
The right panel plots $R(\lambda)$ against $\log_{10}\lambda$. The ratio is
larger than one throughout the range considered, so the Cauchy scale parameter
$\gamma$ is asymptotically more precisely estimated than the Gaussian scale
parameter $\sigma$. This may seem counterintuitive because $\gamma$ governs the
heavy-tailed component. The reason is that the Cauchy scale is strongly
identified by the frequency and magnitude of extreme observations, where the
Gaussian component contributes little. By contrast, the Gaussian scale is
identified primarily from the central part of the distribution, where changes in
$\sigma$ and $\gamma$ can partially offset one another. The U-shape of
$R(\lambda)$ shows that the relative precision is especially favorable to
$\gamma$ when the Cauchy scale is either small or large relative to the Gaussian
scale.
\begin{figure}[ht]
\begin{centering}
\includegraphics[width=1\textwidth]{Figures/SERATIO3.jpg}
\par
\end{centering}
\caption{\small The asymptotic standard deviations for $\sigma$ and $\gamma$ are plotted against the ratio $\lambda=\gamma/\sigma$ in the left panel for $(\sigma,\gamma)=(1,\lambda)$. The right panel shows $R(\lambda)$, the ratio of the asymptotic standard deviation of $\hat\sigma$ to that of $\hat\gamma$. Values above one imply that the Cauchy scale $\gamma$ is estimated more precisely than the Gaussian scale $\sigma$.}
\label{fig:seratio}
\end{figure}
\subsection{Finite-Sample Behavior of the Voigt MLE}
We examine the finite-sample behavior of the Voigt maximum likelihood estimator through a Monte Carlo study. The purpose is to assess how accurately the asymptotic distribution in Theorem \ref{thm:MLE-consistent-asN} approximates the sampling distribution of the MLE. We generate $N=100,000$ independent datasets for sample sizes ranging from $n=100$ to $n=10,000$. For each replication, we compute the MLE and compare the empirical distribution of the estimates with the corresponding asymptotic approximation.
The simulations cover three values of the Cauchy-to-Gaussian scale ratio,
$$
\lambda=\frac{\gamma_0}{\sigma_0},
$$
corresponding to $\lambda=0.01$, $\lambda=0.1$, and $\lambda=1$. Thus the designs range from a nearly Gaussian case to a balanced Voigt case in which the Gaussian and Cauchy scales are equal.
Table \ref{tab:simulation} summarizes the results. We report the mean of the parameter estimates, their empirical standard deviation (Std), and the asymptotic standard deviation (aStd) derived from the Fisher information matrix. For each parameter component, we also report the empirical left and right tail probabilities
\begin{equation}
\alpha_{L,j}
=
\frac{1}{N}
\sum_{i=1}^{N}
\left[
\frac{\hat{\theta}_{n,j}^{(i)}-\theta_{0,j}}
{\operatorname{Std}(\hat{\theta}_{n,j})}
<-1.96
\right],
\qquad
\alpha_{R,j}
=
\frac{1}{N}
\sum_{i=1}^{N}
\left[
\frac{\hat{\theta}_{n,j}^{(i)}-\theta_{0,j}}
{\operatorname{Std}(\hat{\theta}_{n,j})}
>1.96
\right],
\label{eq:RLquantiles}
\end{equation}
where $j$ indexes the components of $\theta=(\mu,\sigma,\gamma)^\prime$. Under asymptotic normality, both probabilities should approach $0.025$ as $n\rightarrow\infty$.
The simulation results support the asymptotic approximation. The estimator of $\mu$ is essentially unbiased across all designs and sample sizes, and its empirical standard deviation is close to the asymptotic standard deviation even for small samples. The scale parameters also behave well once the sample size is moderate. The main finite-sample distortions occur in the smaller samples and involve the scale parameters. In the near-Gaussian design, $\lambda=0.01$, the Cauchy scale is close to the lower boundary of the parameter space, so the finite-sample distribution of $\hat\gamma$ is visibly skewed even though the magnitude of its standard deviation is well approximated by the Fisher information. In the balanced design, $\lambda=1$, the estimator of $\sigma$ exhibits some small-sample skewness, reflecting the difficulty of separating the Gaussian core from the Cauchy tails with limited data. These distortions diminish rapidly as $n$ increases.
Overall, the table indicates that the Voigt MLE is well centered and that the Fisher-information approximation provides an accurate description of sampling variability over a wide range of designs. The remaining discrepancies are finite-sample effects associated with weak identification of one of the two scale components, either because $\gamma_0$ is close to zero or because the sample is too small to clearly distinguish the Gaussian and Cauchy contributions.
\begin{table}[p]
\caption{Simulation Results for the Voigt Density MLE}
\begin{centering}
\vspace{0.2cm}
\begin{footnotesize}
\begin{tabularx}{\textwidth}{Yp{-0.5cm}YYYp{-0.5cm}YYYp{-0.5cm}YYYY}
\toprule
\midrule
& & $\mu$ & $\sigma$ & $\gamma$ & & $\mu$ & $\sigma$ & $\gamma$ & & $\mu$ & $\sigma$ & $\gamma$ \\
\midrule
\\[-0.3cm]
\multicolumn{13}{c}{Panel A: $(\mu_0,\sigma_0, \gamma_0)= (1,1,0.01)$} \\
\\[-0.3cm]
& & \multicolumn{3}{c}{$n=100$} & & \multicolumn{3}{c}{$n=200$} & & \multicolumn{3}{c}{$n=400$} \\
\\[-0.3cm]
Mean & & 0.9999 & 0.9894 & 0.0112 & & 1.0001 & 0.9952 & 0.0102 & & 1.0001 & 0.9980 & 0.0100 \\
Std & & 0.1016 & 0.0791 & 0.0239 & & 0.0715 & 0.0552 & 0.0164 & & 0.0509 & 0.0388 & 0.0117 \\
aStd & & 0.1013 & 0.0775 & 0.0235 & & 0.0716 & 0.0548 & 0.0166 & & 0.0507 & 0.0388 & 0.0118 \\
$\alpha_{L}$ & & 0.0250 & 0.0166 & 0.0588 & & 0.0250 & 0.0199 & 0.0570 & & 0.0258 & 0.0223 & 0.0487 \\
$\alpha_{R}$ & & 0.0249 & 0.0353 & 0.0000 & & 0.0253 & 0.0307 & 0.0000 & & 0.0244 & 0.0290 & 0.0000 \\
\\[-0.2cm]
& & \multicolumn{3}{c}{$n=1,000$} & & \multicolumn{3}{c}{$n=4,000$} & & \multicolumn{3}{c}{$n=10,000$} \\
\\[-0.3cm]
Mean & & 1.0001 & 0.9995 & 0.0099 & & 0.9999 & 0.9998 & 0.0100 & & 0.9999 & 0.9999 & 0.0100 \\
Std & & 0.0320 & 0.0247 & 0.0075 & & 0.0160 & 0.0122 & 0.0037 & & 0.0101 & 0.0077 & 0.0024 \\
aStd & & 0.0320 & 0.0245 & 0.0074 & & 0.0160 & 0.0123 & 0.0037 & & 0.0101 & 0.0078 & 0.0024 \\
$\alpha_{L}$ & & 0.0250 & 0.0233 & 0.0384 & & 0.0252 & 0.0241 & 0.0327 & & 0.0246 & 0.0237 & 0.0301 \\
$\alpha_{R}$ & & 0.0241 & 0.0261 & 0.0000 & & 0.0246 & 0.0258 & 0.0146 & & 0.0254 & 0.0258 & 0.0189 \\
\\[-0.3cm]
\hdashline
\\[-0.3cm]
\multicolumn{13}{c}{Panel B: $(\mu_0,\sigma_0, \gamma_0)= (1,1,0.1)$} \\
\\[-0.3cm]
& & \multicolumn{3}{c}{$n=100$} & & \multicolumn{3}{c}{$n=200$} & & \multicolumn{3}{c}{$n=400$} \\
\\[-0.3cm]
Mean & & 1.0002 & 0.9935 & 0.0972 & & 0.9999 & 0.9970 & 0.0986 & & 1.0003 & 0.9983 & 0.0994 \\
Std & & 0.1122 & 0.1156 & 0.0717 & & 0.0793 & 0.0795 & 0.0505 & & 0.0558 & 0.0552 & 0.0351 \\
aStd & & 0.1112 & 0.1090 & 0.0701 & & 0.0786 & 0.0771 & 0.0496 & & 0.0556 & 0.0545 & 0.0350 \\
$\alpha_{L}$ & & 0.0249 & 0.0202 & 0.0376 & & 0.0252 & 0.0227 & 0.0346 & & 0.0250 & 0.0223 & 0.0322 \\
$\alpha_{R}$ & & 0.0254 & 0.0315 & 0.0000 & & 0.0252 & 0.0294 & 0.0282 & & 0.0247 & 0.0287 & 0.0161 \\
\\[-0.2cm]
& & \multicolumn{3}{c}{$n=1,000$} & & \multicolumn{3}{c}{$n=4,000$} & & \multicolumn{3}{c}{$n=10,000$} \\
\\[-0.3cm]
Mean & & 1.0000 & 0.9994 & 0.0998 & & 1.0000 & 0.9999 & 0.0999 & & 1.0001 & 0.9999 & 0.1000 \\
Std & & 0.0351 & 0.0347 & 0.0223 & & 0.0176 & 0.0173 & 0.0111 & & 0.0111 & 0.0109 & 0.0070 \\
aStd & & 0.0352 & 0.0345 & 0.0222 & & 0.0176 & 0.0172 & 0.0111 & & 0.0111 & 0.0109 & 0.0070 \\
$\alpha_{L}$ & & 0.0257 & 0.0235 & 0.0297 & & 0.0249 & 0.0244 & 0.0272 & & 0.0251 & 0.0238 & 0.0272 \\
$\alpha_{R}$ & & 0.0244 & 0.0270 & 0.0194 & & 0.0254 & 0.0260 & 0.0231 & & 0.0253 & 0.0260 & 0.0234 \\
\\[-0.3cm]
\hdashline
\\[-0.3cm]
\multicolumn{13}{c}{Panel C: $(\mu_0,\sigma_0, \gamma_0)= (1,1,1)$} \\
\\[-0.3cm]
& & \multicolumn{3}{c}{$n=100$} & & \multicolumn{3}{c}{$n=200$} & & \multicolumn{3}{c}{$n=400$} \\
\\[-0.3cm]
Mean & & 0.9996 & 0.9422 & 0.9769 & & 0.9992 & 0.9642 & 0.9929 & & 0.9992 & 0.9846 & 0.9966 \\
Std & & 0.2138 & 0.4571 & 0.2576 & & 0.1492 & 0.3250 & 0.1883 & & 0.1046 & 0.2133 & 0.1338 \\
aStd & & 0.2088 & 0.3910 & 0.2653 & & 0.1477 & 0.2765 & 0.1876 & & 0.1044 & 0.1955 & 0.1326 \\
$\alpha_{L}$ & & 0.0249 & 0.0103 & 0.0206 & & 0.0245 & 0.0086 & 0.0263 & & 0.0248 & 0.0118 & 0.0272 \\
$\alpha_{R}$ & & 0.0255 & 0.0910 & 0.0275 & & 0.0258 & 0.0491 & 0.0232 & & 0.0256 & 0.0357 & 0.0228 \\
\\[-0.2cm]
& & \multicolumn{3}{c}{$n=1,000$} & & \multicolumn{3}{c}{$n=4,000$} & & \multicolumn{3}{c}{$n=10,000$} \\
\\[-0.3cm]
Mean & & 1.0003 & 0.9952 & 0.9985 & & 1.0000 & 0.9992 & 0.9994 & & 1.0000 & 0.9993 & 0.9999 \\
Std & & 0.0662 & 0.1265 & 0.0841 & & 0.0330 & 0.0620 & 0.0418 & & 0.0209 & 0.0393 & 0.0266 \\
aStd & & 0.0660 & 0.1236 & 0.0839 & & 0.0330 & 0.0618 & 0.0419 & & 0.0209 & 0.0391 & 0.0265 \\
$\alpha_{L}$ & & 0.0251 & 0.0180 & 0.0261 & & 0.0249 & 0.0218 & 0.0254 & & 0.0244 & 0.0224 & 0.0253 \\
$\alpha_{R}$ & & 0.0253 & 0.0328 & 0.0244 & & 0.0252 & 0.0278 & 0.0244 & & 0.0251 & 0.0277 & 0.0251 \\
\\[-0.3cm]
\midrule
\bottomrule
\end{tabularx}
\end{footnotesize}
\par\end{centering}
{\small Note: The simulation results are based on $N=100,000$ Monte Carlo replications. We report the mean and standard deviation (Std) of the estimated parameters. For comparison, the asymptotic standard deviations (aStd), computed from the Fisher information, are also included. We also report the empirical left and right tail probabilities, $\alpha_{L}$ and $\alpha_{R}$, defined in (\ref{eq:RLquantiles}); the nominal benchmark for each one-sided tail probability is $0.025$.
\label{tab:simulation}}{\small\par}
\end{table}
\section{The Gauss-Cauchy Convolution Filter}\label{sec:GCCfilter}
We now introduce the Gauss-Cauchy Convolution (GCC) filter, a robust recursive method for extracting a latent state when measurement errors contain both Gaussian background noise and Cauchy outliers. The filter is built from the Voigt likelihood derived above and uses its score, or equivalently the conditional mean of the latent Gaussian component, to update the state
recursively.
The standard Kalman filter is optimal in linear Gaussian systems, but it is not robust to heavy-tailed measurement noise: a single large outlier can move
the state estimate in proportion to its magnitude. The GCC filter replaces the linear Kalman update with a nonlinear update implied by the Gauss-Cauchy
convolution. For moderate prediction errors, the update is approximately linear, resembling the Kalman filter. For large prediction errors, however, the
conditional mean redescends toward zero, so the filter automatically discounts observations that are more plausibly attributed to the Cauchy component than
to the latent Gaussian signal.
A key advantage of the GCC filter is analytical tractability. Conditional on a Gaussian prediction density for the latent state, the prediction-error density is Voigt and therefore available in closed form. This contrasts with many robust alternatives, such as filters based on Student-$t$ measurement noise,
where the convolution of a Gaussian prediction density and a Student-$t$ error does not generally yield a closed-form prediction-error density. The
pure Cauchy measurement-error filter of \citet{McCabeGualdoni:2024} is an important exception and is nested as a limiting case of the present framework.
The GCC filter extends this case by allowing measurement errors to contain
both Gaussian and Cauchy components, a distinction that is empirically useful
when ordinary background noise is approximately Gaussian but occasional
observations are contaminated by heavy-tailed outliers.
Finally, the convolution parameters are estimated by quasi maximum likelihood using
the analytical density, score, and Hessian derived above. Thus the filter
inputs are obtained from likelihood-based estimation rather than ad-hoc
tuning.
\subsection{Filtering Setup}
Let $x_t$ denote a latent state variable that follows a Gaussian AR(1) process, and $y_t$ denote the observed variable, measured with additive error:
\begin{alignat*}{2}
x_t &= (1-\phi)\mu+\phi x_{t-1}+\varepsilon_t,
&\qquad \varepsilon_t &\sim \mathcal{N}(0,\tau^2), \\
y_t &= x_t+\eta_t,
&\qquad \eta_t &\sim \mathcal{V}(0,\sigma,\gamma).
\end{alignat*}
The state innovations, Gaussian measurement-error components, and Cauchy
measurement-error components are mutually independent and independent over
time. The initial state is independent of all subsequent innovations.
Thus the latent state is Gaussian, while the measurement
error follows a Voigt distribution. Equivalently, the measurement error contains
a Gaussian component with scale $\sigma$ and a Cauchy component with scale
$\gamma$.
We use standard state-space notation \citep[see, e.g.,][]{Harvey:1989, DurbinKoopman:2012},
with lowercase letters for the state and observation variables to avoid
conflict with the convolution variables used above. Let
$\mathcal{F}_t=\sigma(y_t,y_{t-1},\ldots)$ denote the natural filtration, and
define
$$
x_{t|t-1}
=
\mathbb{E}[x_t|\mathcal{F}_{t-1}],
\quad\text{and}\quad
h_{t|t-1}
=
\operatorname{var}(x_t|\mathcal{F}_{t-1}).
$$
Here $h_{t|t-1}$ denotes the one-step-ahead prediction variance of the latent
state. The linear Gaussian state transition gives the usual prediction
equations:
$$
x_{t|t-1} = (1-\phi)\mu+\phi x_{t-1|t-1}\quad\text{and}\quad
h_{t|t-1} = \phi^2 h_{t-1|t-1}+\tau^2.
$$
The nonlinear part of the GCC filter enters through the update step. Following
the approximation introduced by \citet{Masreliez:1975}, we approximate the
one-step-ahead prediction density of the state by a Gaussian distribution.
\begin{assumption}[Masreliez approximation]\label{ass:Masreliez}
Conditionally on $\mathcal{F}_{t-1}$,
$$
x_t-x_{t|t-1}\sim\mathcal{N}(0,h_{t|t-1}).
$$
\end{assumption}
This approximation was originally introduced by \citet{Masreliez:1975} for
robust radar tracking in the presence of so-called glint noise, a form of
sporadic heavy-tailed measurement contamination. In the present setting, the
approximation is useful for a more specific reason: it turns the prediction
error into a Gaussian convolution. Let
$$
\xi_t=x_t-x_{t|t-1},
\qquad
e_t=\xi_t+\eta_t.
$$
Under Assumption \ref{ass:Masreliez}, conditionally on $\mathcal F_{t-1}$,
$\xi_t$ is Gaussian and independent of $\eta_t$. Proposition
\ref{prop:TweedieGaussian} can therefore be applied conditionally to obtain the
state correction from the score of the one-step-ahead prediction-error density:
$$
\mathbb{E}[\xi_t|e_t,\mathcal F_{t-1}]
=
-h_{t|t-1}\frac{\partial}{\partial e_t}\log f_{t-1}(e_t),
$$
where $f_{t-1}$ denotes the conditional density of $e_t$ given
$\mathcal F_{t-1}$.
The second ingredient is specific to the GCC model. Since
$\eta_t\sim\mathcal V(0,\sigma,\gamma)$, we may write
$\eta_t=Z_t+C_t$, where $Z_t\sim\mathcal N(0,\sigma^2)$ and
$C_t\sim\operatorname{Cauchy}(0,\gamma)$. Hence
$$
e_t=(\xi_t+Z_t)+C_t,
\qquad
\xi_t+Z_t|\mathcal F_{t-1}\sim\mathcal N(0,h_{t|t-1}+\sigma^2),
$$
so the prediction-error density remains Voigt, specifically
$e_t|\mathcal F_{t-1}\sim\mathcal V(0,\delta_t,\gamma)$ with $\delta_t^2=h_{t|t-1}+\sigma^2$. Thus Assumption \ref{ass:Masreliez} and the Gaussian-Cauchy convolution
structure work together: Tweedie's formula gives the moment update, and the
Voigt closure property makes the required prediction-error score available in
closed form.
\begin{theorem}[GCC filter]\label{thm:GCCfilter}
Under Assumption \ref{ass:Masreliez}, conditionally on $\mathcal{F}_{t-1}$, the prediction error
$e_t=y_t-x_{t|t-1}=x_t-x_{t|t-1}+\eta_t$ satisfies
$$
e_t|\mathcal{F}_{t-1}\sim\mathcal{V}(0,\delta_t,\gamma),
\qquad
\delta_t^2=h_{t|t-1}+\sigma^2.
$$
The conditional mean and variance,
$x_{t|t}=\mathbb{E}[x_t\mid\mathcal{F}_t]$ and
$h_{t|t}=\operatorname{var}(x_t\mid\mathcal{F}_t)$, are updated by
$$
x_{t|t} = x_{t|t-1}+h_{t|t-1}\psi_t
\quad\text{and}\quad
h_{t|t} = h_{t|t-1}-h_{t|t-1}^2\psi_t^\prime,
$$
where
$$
\psi_t = \frac{1}{\delta_t^2}
\left( e_t+\gamma\frac{\mathsf{v}_t}{\mathsf{u}_t}
\right)
\quad\text{and}\quad
\psi_t^\prime
= \frac{1}{\delta_t^4} \left[\delta_t^2+ \gamma^2
\left(1+\frac{\mathsf{v}_t^2}{\mathsf{u}_t^2}
\right)
-\sqrt{\frac{2}{\pi}}\frac{\gamma\delta_t}{\mathsf{u}_t}\right],
$$
with $\psi_t^\prime=\partial\psi_t/\partial e_t$,
$\mathsf{u}_t=\mathsf{u}(e_t;0,\delta_t,\gamma)$, and
$\mathsf{v}_t=\mathsf{v}(e_t;0,\delta_t,\gamma)$.
\end{theorem}
The update has the same structure as the Kalman filter, except that the
linear Gaussian score $e_t/\delta_t^2$ is replaced by the Voigt location score
$\psi_t$. For moderate prediction errors, the update is approximately linear
and resembles the Kalman update. For large prediction errors, however, the
Voigt score redescends toward zero, so the filter discounts observations that
are more plausibly attributed to the Cauchy component than to the latent
Gaussian state.
Figure \ref{fig:score_functions} illustrates this mechanism. The left panel
shows how the Voigt score changes with the Cauchy-to-Gaussian ratio
$\lambda=\gamma/\sigma$. Larger values of $\lambda$ generate stronger
redescending behavior. The right panel compares the GCC score with Gaussian,
Huber, and Student-$t$ scores. The Gaussian score is linear, the Huber score is
bounded, and the Student-$t$ score is redescending. The GCC score is also
redescending, but its shape is tied directly to the Gaussian-Cauchy convolution
rather than imposed as a robustification device.
\begin{figure}[t]
\begin{centering}
\includegraphics[width=1\textwidth]{Figures/ScoreFunctionsPSI.pdf}
\par
\end{centering}
\caption{\small Score functions. The left panel shows the GCC location score, $\psi(e)$,
for different values of $\lambda=\gamma/\sigma$. The right panel compares the
GCC score with Gaussian, Huber, and Student-$t$ scores. The GCC score is
approximately linear near the origin and redescends in the tails, implying that
extreme prediction errors receive little weight in the state update.}
\label{fig:score_functions}
\end{figure}
Under the Gaussian prediction approximation in Assumption \ref{ass:Masreliez},
the recursions in Theorem \ref{thm:GCCfilter} coincide with the approximate
conditional mean filter of \citet{Masreliez:1975}. The update is therefore
optimal within the class of filters that represent the one-step-ahead state
density by a Gaussian approximation and update its first two conditional
moments using the prediction-error density.
The next lemma records that the variance recursion remains positive under the Gaussian prediction approximation.
\begin{lemma}[Positivity of the GCC variance recursion]\label{lem:GCCvariance}
Suppose Assumption \ref{ass:Masreliez} holds and $h_{t|t-1}>0$. Then the GCC
variance update satisfies
$$
h_{t|t}=h_{t|t-1}-h_{t|t-1}^2\psi_t^\prime>0.
$$
Moreover,
$h_{t+1|t} = \phi^2h_{t|t}+\tau^2 \geq \tau^2$,
with strict inequality whenever $\phi\neq0$. Consequently,
$\delta_{t+1}^2=h_{t+1|t}+\sigma^2>0$, so the Voigt prediction-error density remains well defined throughout the
recursion.
\end{lemma}
Next, we introduce the smoothing algorithm.
\subsection{GCC Smoothing}
The GCC filter produces filtered estimates $x_{t|t}$ and $h_{t|t}$. For
applications in which a full-sample estimate of the latent state is desired, we
construct smoothed estimates using a fixed-interval recursion (the Rauch-Tung-Striebel smoother) based on the Gaussian state transition and the moment approximation that underlies the GCC
filter.
Starting from $x_{T|T}$ and $h_{T|T}$, the smoothed estimates are obtained
backward by
$$
x_{t|T}
=
x_{t|t}
+
c_t
\left(
x_{t+1|T}-x_{t+1|t}
\right),
\qquad
c_t
=
\frac{\phi h_{t|t}}{h_{t+1|t}}.
$$
for $t=T-1,\ldots,1$. The smoothed variance is
$$
h_{t|T}
=
h_{t|t}
+
c_t^2
\left(
h_{t+1|T}-h_{t+1|t}
\right).
$$
Here $x_{t|t}$, $h_{t|t}$, $x_{t+1|t}$, and $h_{t+1|t}$ are produced by the
forward GCC filter.
\subsection{Filtering of the Cauchy Components}
Given the system model $y_t = x_t + Z_t + C_t$, where $x_t$ is the latent state, $Z_t \sim \mathcal{N}(0, \sigma^2)$ represents Gaussian measurement noise, and $C_t \sim \operatorname{Cauchy}(0, \gamma)$ captures the heavy-tailed measurement, we can characterize the optimal extraction of the Cauchy component as follows.
Let $\xi_t = x_t - x_{t|t-1}$ denote the state prediction error, and define the total Gaussian unobserved component as $G_t = \xi_t + Z_t$. The prediction error can then be written as $e_t = y_t - x_{t|t-1} = G_t + C_t$.
Under the Masreliez approximation, $G_t \mid \mathcal{F}_{t-1} \sim \mathcal{N}(0, \delta_t^2)$, where the total Gaussian scale is $\delta_t^2 = h_{t|t-1} + \sigma^2$. Applying the exact Voigt identity (Corollary \ref{cor:CondExpect}), the expected total Gaussian component given the concurrent information is:
\begin{equation*}
\mathbb{E}[G_t \mid \mathcal{F}_t] = \delta_t^2 \psi(e_t; 0, \delta_t, \gamma) = e_t + \gamma \frac{\mathsf{v}(e_t; 0, \delta_t, \gamma)}{\mathsf{u}(e_t; 0, \delta_t, \gamma)}.
\end{equation*}
The filtered Cauchy component is isolated residually:
\begin{equation}
\mathbb{E}[C_t \mid \mathcal{F}_t] = e_t - \mathbb{E}[G_t \mid \mathcal{F}_t] = -\gamma \frac{\mathsf{v}(e_t; 0, \delta_t, \gamma)}{\mathsf{u}(e_t; 0, \delta_t, \gamma)}. \label{eq:filtered_cauchy}
\end{equation}
The total Gaussian component $G_t$ is subsequently allocated between the state update and the measurement noise update in proportion to their prior variances:
\begin{equation*}
\mathbb{E}[\xi_t \mid \mathcal{F}_t] = h_{t|t-1}\psi(e_t; 0, \delta_t, \gamma), \qquad \mathbb{E}[Z_t \mid \mathcal{F}_t] = \sigma^2\psi(e_t; 0, \delta_t, \gamma).
\end{equation*}
\subsection{QMLE of the Gauss-Cauchy Convolution Filter}
The parameters of the GCC state-space model are estimated by maximizing the
quasi log-likelihood implied by the Masreliez prediction approximation. Let
$\vartheta=(\mu,\sigma,\gamma,\phi,\tau)^\prime$ and define, recursively,
$e_t(\vartheta)=y_t-x_{t|t-1}(\vartheta)$ and
$\delta_t^2(\vartheta)=h_{t|t-1}(\vartheta)+\sigma^2$. The one-step
quasi log-likelihood contribution is
$$
\ell_t(\vartheta)
=
-\tfrac{1}{2}\log(2\pi)
-\log\delta_t(\vartheta)
+
\log\mathsf{u}(e_t(\vartheta);0,\delta_t(\vartheta),\gamma),
$$
and the full-sample objective is
$$
\ell_T(\vartheta)=\sum_{t=1}^T\ell_t(\vartheta).
$$
The QMLE is
$$
\hat\vartheta_T
=
\arg\max_{\vartheta\in\Theta}\ell_T(\vartheta).
$$
The recursion is initialized at the stationary moments,
$x_{1|0}=\mu$ and $h_{1|0}=\tau^2/(1-\phi^2)$.
The objective is a quasi-likelihood because the prediction-error density is
evaluated under the Gaussian prediction approximation in Assumption
\ref{ass:Masreliez}. Conditional on the recursively generated prediction
moments, however, each likelihood contribution is the exact Voigt density of
the implied prediction error. The analytical derivatives in Lemma
\ref{lem:ScoreHessian} and the algebraic closure property in Corollary
\ref{cor:VoigtAlgebraicClosure} make it possible to compute the likelihood and
its derivatives without numerical convolution or finite differences.
\begin{assumption}
\label{assu:CompactGCC}
The parameter vector $\vartheta=(\mu,\sigma,\gamma,\phi,\tau)^\prime$
belongs to a compact parameter space $\Theta$ satisfying
$|\mu|\leq\mu_{\max}$,
$0<\sigma_{\min}\leq\sigma\leq\sigma_{\max}$,
$0<\gamma_{\min}\leq\gamma\leq\gamma_{\max}$,
$|\phi|\leq\phi_{\max}<1$, and
$0<\tau_{\min}\leq\tau\leq\tau_{\max}$.
\end{assumption}
\begin{assumption}
\label{assu:QMLE}
The pseudo-true parameter
$\vartheta_\star=\arg\max_{\vartheta\in\Theta}
\mathbb{E}_{\vartheta_0}[\ell_t(\vartheta)]$
is unique.
\end{assumption}
\begin{conjecture}[GCC: Asymptotic properties of the QMLE]
\label{conj:MLE-consistent-GCC}
Let $y_1,\ldots,y_T$ be generated by the GCC state-space model with
$\vartheta_0\in\Theta$, and suppose Assumptions
\ref{assu:CompactGCC}-\ref{assu:QMLE} hold. Then
$\hat\vartheta_T \xrightarrow{p} \vartheta_\star$.
If, in addition, $\vartheta_\star$ is an interior point of $\Theta$, then
$$
\frac{1}{\sqrt{T}} \sum_{t=1}^T \frac{\partial \ell_t}{\partial \vartheta} \xrightarrow{d} \mathcal{N} \left(0, \mathcal{I}_{\vartheta_\star}\right) \quad \text { and } \quad -\frac{1}{T} \sum_{t=1}^T \frac{\partial^2 \ell_t}{\partial \vartheta \partial \vartheta^{\prime}} \xrightarrow{p} \mathcal{J}_{\vartheta_\star}
$$
such that
$\sqrt{T}\left(\hat{\vartheta}_T-\vartheta_\star\right) \xrightarrow{d} \mathcal{N} \left(0, \mathcal{J}_{\vartheta_\star}^{-1} \mathcal{I}_{\vartheta_\star} \mathcal{J}_{\vartheta_\star}^{-1} \right)$.
\end{conjecture}
Under the idealized model in which the Gaussian prediction approximation is exact (Assumption \ref{ass:Masreliez}), the quasi-likelihood is correctly specified, so $\vartheta_\star=\vartheta_0$, $\mathcal{I}_{\vartheta_0}=\mathcal{J}_{\vartheta_0}$, and $\sqrt{T}\left(\hat{\vartheta}_T-\vartheta_0\right) \xrightarrow{d} \mathcal{N} \left(0,\mathcal{I}_{\vartheta_0}^{-1} \right)$. For the exact GCC state-space model, we treat this as a conjecture and evaluate it numerically below.
\subsection{QMLE Monte Carlo}
\begin{table}[p]
\caption{Finite-Sample Behavior of the GCC QMLE}
\begin{centering}
\vspace{0.2cm}
\begin{footnotesize}
\begin{tabularx}{\textwidth}{Yp{-0.1cm}YYYYYp{-0.1cm}YYYYY}
\toprule
\midrule
& & $\sigma$ & $\gamma$ & $\mu$ & $\phi$ & $\tau$ & & $\sigma$ & $\gamma$ & $\mu$ & $\phi$ & $\tau$ \\
\midrule
\\[-0.3cm]
\multicolumn{13}{c}{Panel A: $(\sigma_0, \gamma_0,\mu_0,\phi_0,\tau_0)= (1,0.1,1,0.95,1)$} \\
\\[-0.3cm]
& & \multicolumn{5}{c}{$T=500$} & & \multicolumn{5}{c}{$T=1,000$} \\
\\[-0.3cm]
Mean & & 0.9889 & 0.0994 & 1.0007 & 0.9389 & 1.0055 & & 0.9945 & 0.0998 & 0.9987 & 0.9446 & 1.0027 \\
Std & & 0.0986 & 0.0376 & 0.8401 & 0.0190 & 0.0870 & & 0.0684 & 0.0264 & 0.6151 & 0.0123 & 0.0611 \\
aStd & & 0.0947 & 0.0368 & 0.8927 & 0.0158 & {0.0876} & & 0.0670 & 0.0260 & 0.6313 & {0.0112} & {0.0620} \\
$\alpha_{L}$ & & 0.0128 & 0.0326 & 0.0244 & 0.0000 & 0.0298 & & 0.0167 & 0.0312 & 0.0248 & 0.0005 & 0.0282 \\
$\alpha_{R}$ & & 0.0383 & 0.0149 & 0.0251 & 0.0928 & 0.0208 & & 0.0349 & 0.0181 & 0.0247 & 0.0751 & 0.0213 \\
\\[-0.2cm]
& & \multicolumn{5}{c}{$T=2,000$} & & \multicolumn{5}{c}{$T=4,000$} \\
\\[-0.3cm]
Mean & & 0.9973 & 0.0998 & 0.9997 & 0.9474 & 1.0013 & & 0.9985 & 0.1000 & 1.0001 & 0.9487 & 1.0006 \\
Std & & 0.0478 & 0.0187 & 0.4416 & 0.0083 & 0.0432 & & 0.0335 & 0.0132 & 0.3133 & 0.0057 & 0.0303 \\
aStd & & 0.0474 & 0.0184 & 0.4464 & 0.0079 & {0.0438} & & 0.0335 & 0.0130 & 0.3156 & {0.0056} & {0.0310} \\
$\alpha_{L}$ & & 0.0189 & 0.0293 & 0.0251 & 0.0028 & 0.0267 & & 0.0210 & 0.0277 & 0.0256 & 0.0069 & 0.0262 \\
$\alpha_{R}$ & & 0.0322 & 0.0205 & 0.0245 & 0.0607 & 0.0229 & & 0.0303 & 0.0213 & 0.0250 & 0.0499 & 0.0240 \\
\\[-0.2cm]
& & \multicolumn{5}{c}{$T=8,000$} & & \multicolumn{5}{c}{$T=16,000$} \\
\\[-0.3cm]
Mean & & 0.9991 & 0.1001 & 1.0000 & 0.9494 & 1.0003 & & 0.9993 & 0.1001 & 0.9993 & 0.9497 & 1.0001 \\
Std & & 0.0236 & 0.0093 & 0.2216 & 0.0040 & 0.0214 & & 0.0167 & 0.0066 & 0.1572 & 0.0028 & 0.0151 \\
aStd & & 0.0237 & 0.0092 & 0.2232 & 0.0039 & {0.0219} & & 0.0167 & 0.0065 & 0.1578 & {0.0028} & {0.0155} \\
$\alpha_{L}$ & & 0.0214 & 0.0273 & 0.0254 & 0.0114 & 0.0259 & & 0.0213 & 0.0277 & 0.0250 & 0.0152 & 0.0255 \\
$\alpha_{R}$ & & 0.0281 & 0.0219 & 0.0244 & 0.0425 & 0.0243 & & 0.0281 & 0.0229 & 0.0252 & 0.0377 & 0.0247 \\
\\[-0.3cm]
\hdashline
\\[-0.3cm]
\multicolumn{13}{c}{Panel B: $(\sigma_0, \gamma_0,\mu_0,\phi_0,\tau_0)= (1,1,1,0.95,1)$} \\
\\[-0.3cm]
& & \multicolumn{5}{c}{$T=500$} & & \multicolumn{5}{c}{$T=1,000$} \\
\\[-0.3cm]
Mean & & 0.9295 & 0.9938 & 1.0013 & 0.9377 & {1.0152} & & 0.9620 & 0.9988 & 0.9989 & {0.9441} & {1.0084} \\
Std & & 0.3602 & 0.1358 & 0.8455 & 0.0223 & {0.1308} & & 0.2462 & 0.0970 & 0.6170 & {0.0141} & {0.0912} \\
aStd & & 0.3067 & 0.1378 & 0.8929 & 0.0177 & {0.1301} & & 0.2169 & 0.0974 & 0.6314 & {0.0125} & {0.0920} \\
$\alpha_{L}$ & & 0.0024 & 0.0233 & 0.0246 & 0.0000 & {0.0360} & & 0.0043 & 0.0263 & 0.0250 & {0.0003} & {0.0333} \\
$\alpha_{R}$ & & 0.0692 & 0.0258 & 0.0252 & 0.0882 & {0.0164} & & 0.0472 & 0.0239 & 0.0248 & {0.0743} & {0.0182} \\
\\[-0.2cm]
& & \multicolumn{5}{c}{$T=2,000$} & & \multicolumn{5}{c}{$T=4,000$} \\
\\[-0.3cm]
Mean & & 0.9801 & 0.9998 & 0.9993 & 0.9471 & {1.0049} & & 0.9878 & 1.0004 & 1.0001 & {0.9485} & {1.0031} \\
Std & & 0.1617 & 0.0688 & 0.4410 & 0.0094 & {0.0637} & & 0.1108 & 0.0488 & 0.3113 & {0.0065} & {0.0448} \\
aStd & & 0.1534 & 0.0689 & 0.4465 & 0.0088 & {0.0651} & & 0.1084 & 0.0487 & 0.3157 & {0.0062} & {0.0460} \\
$\alpha_{L}$ & & 0.0094 & 0.0264 & 0.0255 & 0.0024 & {0.0312} & & 0.0132 & 0.0273 & 0.0262 & {0.0064} & {0.0312} \\
$\alpha_{R}$ & & 0.0423 & 0.0234 & 0.0253 & 0.0620 & {0.0199} & & 0.0396 & 0.0232 & 0.0260 & {0.0516} & {0.0213} \\
\\[-0.2cm]
& & \multicolumn{5}{c}{$T=8,000$} & & \multicolumn{5}{c}{$T=16,000$} \\
\\[-0.3cm]
Mean & & 0.9913 & 1.0006 & 0.9998 & 0.9493 & {1.0023} & & 0.9928 & 1.0005 & 0.9997 & {0.9496} & {1.0020} \\
Std & & 0.0774 & 0.0343 & 0.2178 & 0.0045 & {0.0315} & & 0.0540 & 0.0242 & 0.1517 & {0.0031} & {0.0221} \\
aStd & & 0.0767 & 0.0344 & 0.2232 & 0.0044 & {0.0325} & & 0.0542 & 0.0244 & 0.1579 & {0.0031} & {0.0230} \\
$\alpha_{L}$ & & 0.0141 & 0.0281 & 0.0268 & 0.0108 & {0.0320} & & 0.0156 & 0.0280 & 0.0275 & {0.0146} & {0.0326} \\
$\alpha_{R}$ & & 0.0375 & 0.0237 & 0.0265 & 0.0431 & {0.0203} & & 0.0382 & 0.0244 & 0.0275 & {0.0394} & {0.0202} \\
\\[-0.3cm]
\midrule
\bottomrule
\end{tabularx}
\end{footnotesize}
\par\end{centering}
{\small Note: The table reports Monte Carlo results for the QMLE of the GCC filtering model, based on $N=100,000$ replications. We report the mean and standard deviation (Std) of the parameter estimates. The asymptotic standard deviations (aStd) are computed from the Fisher information using a long simulated sample with $T=10,000,000$ observations. The empirical left and right tail probabilities, $\alpha_L$ and $\alpha_R$, are defined as in (\ref{eq:RLquantiles}); the nominal benchmark for each one-sided probability is 0.025.
\label{tab:GCC_QMLE_simulation}}{\small\par}
\end{table}
Table \ref{tab:GCC_QMLE_simulation} reports Monte Carlo evidence for the
finite-sample behavior of the GCC QMLE. The simulations cover two measurement-error designs: a small Cauchy component, $\gamma_0/\sigma_0=0.1$, and a balanced design, $\gamma_0/\sigma_0=1$. In both cases the estimators are well centered, and the empirical standard deviations are close to the inverse-information standard deviations even for moderate sample sizes.
The small-Cauchy design in Panel A is especially regular. The estimates of
$\sigma$, $\gamma$, $\mu$, and the state innovation scale are nearly unbiased throughout the table, and the empirical standard deviations track the asymptotic standard deviations closely. The main finite-sample distortion is in the persistence parameter $\phi$, which is biased downward in small samples and has right-tail rejection probabilities above the nominal benchmark. This distortion shrinks steadily as $T$ increases, although at sample sizes typical of financial applications the remaining bias in $\hat\phi$ can still imply some understatement of persistence. This is most naturally interpreted as a finite-sample dynamic-parameter issue rather than a Masreliez approximation effect.
Panel B shows the more challenging case in which the Gaussian and Cauchy
measurement-error scales are equal. The estimator of $\gamma$ remains well
centered even in small samples, whereas $\hat\sigma$ displays a noticeable
downward bias for small $T$. This reflects the difficulty of separating the
Gaussian core from the Cauchy component when the two scales are of comparable
magnitude and the sample is short. The bias decreases monotonically with $T$,
and by $T=8{,}000$ or $T=16{,}000$ the empirical standard deviations are very
close to their asymptotic counterparts.
The tail probabilities provide a useful diagnostic for the normal approximation.
For $\mu$, the empirical left and right probabilities are close to the nominal
$0.025$ benchmark across both panels, even for small samples. The scale and
persistence parameters show more finite-sample skewness, especially $\sigma$ in
Panel B and $\phi$ in both panels. These deviations diminish with the sample
size, consistent with the regular inverse-information approximation becoming
accurate away from boundary cases.
\begin{comment}
Overall, the results support the use of the recursively evaluated quasi-likelihood
for estimation of the GCC filter. The QMLE is well centered, the inverse-information
standard errors provide a good approximation to sampling variability, and the
remaining discrepancies are concentrated in the small-sample behavior of the
scale decomposition and the persistence parameter.
\end{comment}
The Monte Carlo evidence is consistent with the conjecture that the pseudo-true parameter, $\vartheta_{\star}$, coincides with the true data-generating parameter, $\vartheta_{0}$. Across the designs considered, the QMLE is well centered, and the inverse-information standard errors provide a good approximation to the sampling variability. This suggests that the information loss induced by the Masreliez approximation is small in these settings. The remaining discrepancies are concentrated in the small-sample behavior of the scale decomposition and the persistence parameter.
\begin{comment}
{\color{red} The fact that the QMLE is well centered in finite samples provides strong numerical evidence that $\hat{\vartheta}_T \xrightarrow{p} \vartheta_0$ (i.e., the pseudo-true parameter $\vartheta_{\star}$ effectively coincides with the true data-generating parameter $\vartheta_0$ ). That the estimator is consistent for the true parameters and the inverse-information standard errors accurately capture the sampling variability collectively demonstrates that the information loss from the Masreliez approximation is negligible. The remaining discrepancies are concentrated in the small-sample behavior of the scale decomposition and the persistence parameter.}
\end{comment}
\section{Assessing the Masreliez Approximation in the GCC Filter}
\label{sec:Simulations}
The GCC filter rests on two linked implications of Assumption
\ref{ass:Masreliez}. First, the Gaussian approximation to the predictive state
density makes the update a Tweedie moment correction for a Gaussian convolution.
Second, conditional on this approximation, the GCC measurement-error structure
keeps the prediction-error density Voigt, so the required score is available in
closed form. The second step is exact given the approximation; the simulations
below assess the numerical cost of the first step.
The purpose of this section is therefore not to evaluate finite-sample maximum
likelihood performance. Instead, we hold the model parameters fixed at their
true values and assess the numerical consequences of Assumption
\ref{ass:Masreliez} for the filtering recursion itself. The central question is
how much is lost by replacing the exact one-step-ahead predictive distribution
of the latent state by the Gaussian approximation used by the GCC filter.
The benchmark is an exact filter computed by numerical density propagation, in
the spirit of the non-Gaussian filtering and smoothing approach of
\citet{Kitagawa:1987}. At each date, Bayes' rule is used to update the full
predictive density of the latent state, and the resulting filtering density is
propagated forward through the Gaussian state transition. This benchmark does
not impose Gaussianity on the predictive state density. It is therefore useful
for assessing the accuracy of the Masreliez approximation, although it is not
the operational filter we propose. Details on the numerical benchmark and the
diagnostic definitions are provided in Appendix \ref{app:simulation-details}.
We compare the exact benchmark with two Gaussian approximations. The first is a
moment-matched Gaussian approximation that uses the exact predictive mean and
variance. This isolates the shape error from Gaussianizing the exact predictive
density. The second is the operational GCC approximation, which uses the
predictive mean and variance generated by the GCC recursion. This measures the
full approximation error of the implemented filter.
We evaluate the approximation in two ways. First, we compute Kullback-Leibler
discrepancies for both the latent-state predictive density and the implied
observation predictive density. Second, we compare the one-step correction
generated by the exact benchmark with the correction generated by the GCC
filter, since this correction is the object that drives the filtered state.
Because the model is equivariant to common changes in scale, the simulation
design is indexed by the dimensionless ratios $\lambda=\gamma/\sigma$ and
$\tau/\sigma$, together with the persistence parameter $\phi$. The grid includes
$\lambda=0.10$, which is close to the value estimated in the empirical
application.
\begin{table}[tbp]
\centering
\small
\caption{Kullback-Leibler diagnostics for the Gaussian prediction approximation}
\label{tab:gcc_density_diagnostics_main}
\vspace{0.2cm}
\begin{tabularx}{\textwidth}{@{} YYYYYY @{}}
\toprule
\midrule
$\lambda$ & $\overline{\mathrm{KL}}^{x,\operatorname{shape}}$ & $\overline{\mathrm{KL}}^{x,\operatorname{op}}$ & $\overline{\mathrm{KL}}^{y,\operatorname{shape}}$ & $\overline{\mathrm{KL}}^{y,\operatorname{op}}$ & $\max\mathrm{KL}^{x,\operatorname{op}}$ \\
\midrule
\\[-0.2cm]
$0.00$ & $<10^{-15}$ & $<10^{-15}$ & $<10^{-15}$ & $<10^{-15}$ & $<10^{-15}$ \\
$0.01$ & $7.05\times 10^{-5}$ & $9.10\times 10^{-5}$ & $1.34\times 10^{-5}$ & $2.17\times 10^{-5}$ & $2.78\times 10^{-2}$ \\
$0.05$ & $1.78\times 10^{-4}$ & $2.33\times 10^{-4}$ & $2.82\times 10^{-5}$ & $4.76\times 10^{-5}$ & $2.13\times 10^{-2}$ \\
$0.10$ & $2.94\times 10^{-4}$ & $4.11\times 10^{-4}$ & $3.93\times 10^{-5}$ & $8.38\times 10^{-5}$ & $3.50\times 10^{-2}$ \\
$0.50$ & $7.30\times 10^{-4}$ & $1.03\times 10^{-3}$ & $4.57\times 10^{-5}$ & $1.31\times 10^{-4}$ & $6.47\times 10^{-2}$ \\
$1.00$ & $8.44\times 10^{-4}$ & $1.24\times 10^{-3}$ & $2.88\times 10^{-5}$ & $1.06\times 10^{-4}$ & $4.11\times 10^{-2}$ \\
\\[-0.2cm]
\\[-0.5cm]
\midrule
\bottomrule
\end{tabularx}
\smallskip\par
\begin{minipage}{\textwidth}
\small
\textit{Note:} The table reports density-level Kullback-Leibler diagnostics for the Masreliez-GCC prediction approximation. For each value of $\lambda=\gamma/\sigma$, entries are averaged over the nine designs generated by $\phi\in\{0.90,0.97,0.99\}$ and $\tau/\sigma\in\{0.25,0.50,1.00\}$. The superscript $x$ refers to the predictive density of the latent state, while $y$ refers to the predictive density of the observation. The ``shape'' columns compare the exact predictive density with its moment-matched Gaussian approximation, whereas the ``op'' columns compare the exact predictive density with the operational Gaussian approximation used by the GCC filter. The final column reports the maximum value of $\mathrm{KL}^{x,\operatorname{op}}$ over all dates and designs within each row.
\end{minipage}
\end{table}
\begin{table}[tbp]
\centering
\small
\caption{Correction distortion by Cauchy-to-Gaussian measurement-noise ratio}
\label{tab:gcc_correction_diagnostics_main}
\vspace{0.2cm}
\begin{tabularx}{\textwidth}{@{}YYYYY@{}}
\toprule
\midrule
$\lambda$
& $\operatorname{MAE}^{\operatorname{shape}}$
& $\operatorname{MAE}^{\operatorname{op}}$
& $\operatorname{RMSE}^{\operatorname{op}}$
& $q_{.95}|D^{\operatorname{op}}|$ \\
\midrule
\\[-0.2cm]
$0.00$ & $<10^{-15}$ & $<10^{-15}$ & $<10^{-15}$ & $<10^{-15}$ \\
$0.01$ & $6.53\times 10^{-4}$ & $1.10\times 10^{-3}$ & $4.72\times 10^{-3}$ & $3.83\times 10^{-3}$ \\
$0.05$ & $1.84\times 10^{-3}$ & $2.80\times 10^{-3}$ & $7.71\times 10^{-3}$ & $1.06\times 10^{-2}$ \\
$0.10$ & $3.41\times 10^{-3}$ & $4.94\times 10^{-3}$ & $1.28\times 10^{-2}$ & $1.96\times 10^{-2}$ \\
$0.50$ & $8.24\times 10^{-3}$ & $1.14\times 10^{-2}$ & $2.10\times 10^{-2}$ & $3.81\times 10^{-2}$ \\
$1.00$ & $9.18\times 10^{-3}$ & $1.32\times 10^{-2}$ & $2.17\times 10^{-2}$ & $4.34\times 10^{-2}$ \\
\\[-0.2cm]
\\[-0.5cm]
\midrule
\bottomrule
\end{tabularx}
\smallskip\par
\begin{minipage}{\textwidth}
\small
\textit{Note:} The table reports correction-level diagnostics for the Masreliez-GCC approximation. For each value of $\lambda=\gamma/\sigma$, entries are averaged over the nine designs generated by $\phi\in\{0.90,0.97,0.99\}$ and $\tau/\sigma\in\{0.25,0.50,1.00\}$. The distortion is $D_t(y_t)=\Delta_t^\ast(y_t)-\Delta_t(y_t)$, where $\Delta_t^\ast(y_t)$ is the exact correction and $\Delta_t(y_t)$ is the operational GCC correction. The shape diagnostic replaces the exact predictive density by its moment-matched Gaussian approximation.
\end{minipage}
\end{table}
The density-level diagnostics in Table \ref{tab:gcc_density_diagnostics_main}
show that the Gaussian prediction approximation is highly accurate across the
simulation designs. For empirically relevant values $\lambda\leq0.10$, the
average KL discrepancies for the latent-state density are on the order of
$10^{-4}$ or smaller, and the discrepancies for the observation density are an
order of magnitude smaller. This attenuation reflects the smoothing effect of
convolution with the Voigt measurement-error density.
Table \ref{tab:gcc_correction_diagnostics_main} shows that these density-level
differences have little effect on the filter update. For $\lambda\leq0.10$, the
mean absolute operational distortion is between $10^{-3}$ and $5\times10^{-3}$,
with root mean square errors below $1.5\times10^{-2}$. The operational
diagnostics are only modestly larger than the shape diagnostics, indicating that
the GCC recursion does not materially amplify the approximation error.
The largest discrepancies occur in stress designs with $\lambda=0.50$ or
$\lambda=1.00$, but even there the distortions remain moderate. Overall, the
evidence suggests that violations of Assumption \ref{ass:Masreliez} have limited
practical consequences for the GCC filter, especially in the empirically
relevant region where the Cauchy component is small relative to the Gaussian
measurement-error component.
\section{Empirical Application}\label{sec:Empirical}
We now apply the GCC filter to realized volatility data. The objective is to
evaluate whether the Gauss-Cauchy measurement-error specification provides a
useful empirical decomposition of realized volatility into a persistent latent
component and transitory measurement noise. We compare the GCC filter with
several alternative filters that use the same Gaussian AR(1) state equation but
differ in their measurement-error distributions.
\subsection{Data}
We apply the GCC filter to the daily logarithm of realized volatility for the
Technology Select Sector SPDR Fund (XLK). The sample runs from December 22,
1998 to November 7, 2025, and the realized-volatility series is obtained from
the Dacheng Xiu Risk Lab. Realized measures are constructed from high-frequency
intraday returns and are therefore informative about daily variation in
volatility, but they may also contain extreme observations associated with
market stress, liquidity disruptions, and market microstructure effects.
We work with log realized volatility. This transformation reduces the
right-skewness of realized volatility and makes a symmetric measurement-error
specification more plausible. In the XLK sample, the skewness of the log
realized-volatility series is 0.084.
\subsection{Competing filters}
We compare the GCC filter with several filters that retain the same Gaussian
AR(1) state equation but impose different specifications on the measurement
error $\eta_t$. The comparison is designed to separate two issues. The first is
whether the measurement-error distribution is sufficiently flexible to account
for outliers in realized volatility. The second is whether the implied
prediction-error density is available in closed form after convolution with the
Gaussian state prediction error.
Table \ref{tab:filter_comparison} summarizes the competing specifications. The
Gaussian filter corresponds to the Kalman benchmark. The Cauchy filter is the
opposite limiting case in which all measurement error is heavy-tailed. The
Normal-Laplace filter provides a light-tailed robust alternative with a
closed-form Gaussian convolution. The Student-$t$ and Huber filters are included
as common robust benchmarks, but their prediction-error densities are not closed
under convolution with Gaussian state uncertainty. Implementations of these
filters therefore require either numerical integration or an approximation to
the prediction-error density.
\begin{table}[tbp]
\centering
\small
\caption{Measurement-error specifications used in the empirical comparison}
\label{tab:filter_comparison}
\vspace{0.2cm}
\begin{tabularx}{\textwidth}{@{} l >{\hsize=0.9\hsize}Y >{\hsize=1.1\hsize}Y c @{}}
\toprule
\midrule
Filter & Measurement error $\eta_t$ & Prediction-error density for $e_t$ & Closed form \\
\midrule
\\[-0.2cm]
Gaussian & $\mathcal{N}(0,\sigma^2)$
& Gaussian with variance $h_{t|t-1}+\sigma^2$ & Yes \\
Cauchy & $\operatorname{Cauchy}(0,\gamma)$
& Voigt convolution of $\mathcal{N}(0,h_{t|t-1})$ and Cauchy error & Yes \\
GCC & $\mathcal{V}(0,\sigma,\gamma)$
& Voigt with Gaussian scale $\delta_t^2=h_{t|t-1}+\sigma^2$ and Cauchy scale $\gamma$ & Yes \\
Normal-Laplace & $\mathcal{N}(0,\sigma^2)+\operatorname{Laplace}(0,\gamma)$
& Normal-Laplace with Gaussian scale $h_{t|t-1}+\sigma^2$ & Yes \\
Student-$t$ & $t_\nu(0,\sigma)$
& Gaussian-$t$ convolution; not generally a Student-$t$ density & No \\
Huber & Huber density with Gaussian core and exponential tails
& Gaussian-Huber convolution; no elementary closed form & No \\
\\[-0.2cm]
\\[-0.5cm]
\midrule
\bottomrule
\end{tabularx}
\smallskip\par
\begin{minipage}{\textwidth}
\small
\textit{Note:} All filters use the same Gaussian AR(1) state equation. The
distinction is the measurement-error distribution and whether convolution with
the Gaussian prediction error preserves an analytically tractable
prediction-error density. Additional details on the parameterizations of the
Normal-Laplace, Huber, and Student-$t$ filters are provided in
Appendix~\ref{app:filter_details}.
\end{minipage}
\end{table}
\subsection{Estimation results}
Table \ref{tab:DFest-block} reports estimation results for the competing
filters. The GCC specification attains the largest value of the implemented
prediction-error criterion, $\ell=-1,306$, substantially improving upon the
Gaussian benchmark, $\ell=-2,465$. This large difference indicates that a
purely Gaussian measurement-error specification is too restrictive for log
realized volatility.
Allowing for heavy tails improves the criterion across all alternatives. The Student-$t$ filter increases the implemented criterion to $\ell=-1,349$, with an estimated degrees-of-freedom parameter $\nu\approx 5.25$, while the pure Cauchy specification yields $\ell=-1,536$. The comparison with the Student-$t$ and Huber filters should be interpreted as a comparison of the operational filtering rules, since their reported values are pseudo-log likelihoods based on same-family approximations to the prediction-error density. With this qualification, the empirical ranking favors the GCC filter.
The Normal-Laplace estimate is also informative: its Gaussian scale is essentially zero, so the fitted specification is close to a pure-Laplace measurement-error model. Thus, within the Normal-Laplace family, the data do not support a separately identified Gaussian measurement-error component once the Laplace component is included.
The GCC estimates imply a measurement-error decomposition with a dominant
Gaussian component and a small but non-negligible Cauchy component. The
estimated parameters are $\sigma=0.1817$ and $\gamma=0.0199$, corresponding to
$\lambda=\gamma/\sigma\approx 0.11$. This estimate lies in the range for which the simulation results in Section \ref{sec:Simulations} show that the Gaussian prediction approximation is essentially exact. Thus, the improvement is not driven by replacing Gaussian noise with a purely heavy-tailed specification, but by allowing for occasional large deviations while retaining a Gaussian core. In contrast, the Gaussian filter absorbs these extreme observations by inflating the measurement variance, $\sigma=0.2980$, which leads to a noisier extraction of the latent state.
Beyond empirical fit, the GCC specification also has an analytical advantage.
As indicated in Table \ref{tab:DFest-block}, the prediction-error density is available in closed form for the GCC, Gaussian, Cauchy, and Normal-Laplace models under the Gaussian prediction approximation. This allows the likelihood contribution and score to be evaluated without numerical convolution. By contrast, the Student-$t$ and Huber filters are not closed under convolution with Gaussian state uncertainty, so their implementations rely on approximations to the prediction-error density.
Finally, all models imply highly persistent latent volatility dynamics, with estimates of $\phi$ close to $0.97$. The exception is the pure Cauchy
specification, which yields a somewhat lower persistence estimate
($\phi\approx 0.936$), reflecting its greater reliance on the measurement-error component to absorb large fluctuations.
\begin{table}[htbp]
\caption{Estimation Results under Different Filtering Methods}
\begin{centering}
\vspace{0.2cm}
\begin{small}
\begin{tabularx}{\textwidth}{p{1.95cm}YYYYYYY}
\toprule
\midrule
\multicolumn{1}{c}{Filter}
& \begin{tabular}{c}GCC\\[-0.5mm](Voigt)\end{tabular}
& \begin{tabular}{c}Kalman\\[-0.5mm](Gaussian)\end{tabular}
& Cauchy
& \begin{tabular}{c}Normal-\\[-0.5mm]Laplace\end{tabular}
& Student-$t$
& Huber \\
\midrule
\\[-0.2cm]
\multicolumn{1}{c}{$\sigma$} & 0.1817 & 0.2980 & & 0.0002 & 0.1696 & 0.1755 \\
& {\it (0.0053)} & {\it (0.0109)} & & {\it (0.0001)} & {\it (0.0044)} & {\it (0.0135)} \\
\multicolumn{1}{c}{$\gamma$} & 0.0199 & & 0.0553 & 0.1786 & & \\
& {\it (0.0021)} & & {\it (0.0037)} & {\it (0.0042)} & & \\
\multicolumn{1}{c}{$\nu$} & & & & & 5.2545 & \\
& & & & & {\it (0.3289)} & \\
\multicolumn{1}{c}{$k$ } & & & & & & 1.3231 \\
& & & & & & {\it (0.1931)} \\
\multicolumn{1}{c}{$\mu$ } & -1.9420 & -1.9492 & -1.9529 & -1.9600 & -1.9667 & -1.9386 \\
& {\it (0.0502)} & {\it (0.0535)} & {\it (0.0351)} & {\it (0.0490)} & {\it (0.0461)} & {\it (0.0622)} \\
\multicolumn{1}{c}{$\phi$} & 0.9716 & 0.9768 & 0.9361 & 0.9735 & 0.9699 & 0.9681 \\
& {\it (0.0040)} & {\it (0.0037)} & {\it (0.0070)} & {\it (0.0034)} & {\it (0.0039)} & {\it (0.0058)} \\
\multicolumn{1}{c}{$\tau$} & 0.1138 & 0.1010 & 0.1813 & 0.1087 & 0.1145 & 0.1183 \\
& {\it (0.0055)} & {\it (0.0086)} & {\it (0.0072)} & {\it (0.0042)} & {\it (0.0053)} & {\it (0.0052)} \\
& & & & & & \\
\multicolumn{1}{c}{$\operatorname{var}(x_t)$} & 0.2316 & 0.2236 & 0.2659 & 0.2257 & 0.2209 & 0.2239 \\
\\[-0.2cm]
\multicolumn{1}{c}{$\ell$} & \textbf{-1,306} & -2,465 & -1,536 & -1,474 & -1,349 & -1,444 \\
& & & & & & \\
Closed-form & Yes & Yes & Yes & Yes & No & No \\
\\[-0.2cm]
\\[-0.5cm]
\midrule
\bottomrule
\end{tabularx}
\end{small}
\par\end{centering}
{\small \textit{Note:} Estimation results under different filtering methods for the logarithm of daily realized volatility for the Technology Select Sector SPDR Fund (XLK). Sample: December 22, 1998 to November 7, 2025. Column headers indicate the filter name, with the measurement-error distribution in parentheses when it differs from the filter name. The parameter $k$ denotes the standardized threshold in the Huber density. $\operatorname{var}(x_t)=\tau^2/(1-\phi^2)$ is the unconditional variance of the latent state. We also report the maximized full-sample prediction-error criterion, $\ell$, and indicate the largest value in bold. For GCC, Kalman, Cauchy, and Normal-Laplace, this criterion is the log-likelihood based on the closed-form prediction-error density under the Gaussian prediction approximation. For Student-$t$ and Huber, it is a pseudo-log likelihood based on the standard same-family approximation to the prediction-error density. The last row indicates whether the filtering method has a closed-form expression for the prediction-error density. Sandwich-type robust standard errors are provided in parentheses. \label{tab:DFest-block}}
\end{table}
\subsection{Filtered volatility and residual decomposition}
Figure \ref{fig:logrv} plots the observed log realized-volatility series and
the filtered GCC estimate of the latent volatility component. The observed
series contains several large spikes associated with periods of market stress,
including the bursting of the Dot-Com bubble, the period following September
11, 2001, the Flash Crash on May 6, 2010, and the market disruption on August
24, 2015. The sample also contains unusually low realized-volatility
measurements in the early part of the sample, which are plausibly related to
market microstructure effects and thin trading.
The filtered state is substantially smoother than the observed series and is
not pulled strongly toward isolated extreme observations. This behavior is a
direct consequence of the redescending location score of the Voigt
prediction-error density, for which $\psi_t\to0$ as $|e_t|\to\infty$. Large
prediction errors are therefore downweighted in the state update and are instead
attributed primarily to the heavy-tailed component of the measurement error.
The right panel of Figure \ref{fig:logrv}, labeled ``Measurement errors,''
illustrates the filtered measurement-error decomposition. The Gaussian
component captures regular, small-to-moderate measurement variation around the
latent state, while the Cauchy component absorbs the most extreme deviations.
Thus, the GCC filter preserves sensitivity to ordinary daily variation in
realized volatility while limiting the influence of observations that are more
plausibly interpreted as transitory measurement disturbances.
\begin{figure}[t]
\begin{centering}
\includegraphics[width=1\textwidth]{Figures/GCC2.jpg}
\par
\end{centering}
\caption{\small GCC filtering results for the daily logarithm of realized
volatility of the Technology Select Sector SPDR Fund (XLK), December 22, 1998
to November 7, 2025. The left panel shows the observed series $y_t$ together
with the filtered state estimate. The right panel presents the filtered
\ Gaussian and Cauchy measurement-error components by using (\ref{eq:filtered_cauchy}).}
\label{fig:logrv}
\end{figure}
\section{Conclusion}\label{sec:Conclusion}
This paper develops likelihood-based inference and filtering methods for the Gauss-Cauchy convolution model. The resulting Voigt density has long been known in spectroscopy, but its statistical use has been limited by the perception that likelihood evaluation, differentiation, and filtering are computationally burdensome. We show that this perception is misleading.
The key observation is that the scaled complementary error function provides an algebraically tractable representation of the Voigt density. Its differential properties yield closed-form expressions for the likelihood, score, Hessian, and conditional moments. As a result, maximum likelihood estimation for the Gauss-Cauchy convolution can be implemented without numerical convolution, quadrature, finite-difference derivatives, or pseudo-Voigt approximations.
The same structure also yields a natural robust filtering method. In the GCC state-space model, the prediction error has a Voigt density under the Masreliez Gaussian prediction approximation. The state update is therefore driven by the location score of this density. This score is approximately
linear near the center but redescends in the tails, so extreme observations are discounted rather than treated as persistent movements in the latent state. In this sense, the GCC filter provides a probabilistic alternative to ad hoc robust filtering rules: robustness arises from the conditional expectation implied by the convolution model itself.
We also assess the numerical consequences of the Masreliez approximation directly. Using an exact benchmark filter based on density propagation, we show that the Gaussian prediction approximation is highly accurate over a broad
range of parameter configurations. The approximation error is small both at the level of predictive densities and, more importantly, at the level of the one-step correction that drives the filter. In the empirically relevant region
where the estimated Cauchy-to-Gaussian ratio is small, the approximation error is negligible.
The empirical application to log realized volatility for XLK illustrates the practical value of the approach. The GCC filter improves substantially on the Gaussian benchmark and outperforms common robust alternatives. The estimates
indicate a dominant Gaussian measurement-error component together with a small Cauchy component that absorbs occasional extreme observations. This decomposition allows the filtered volatility state to remain smooth without
being unduly influenced by transient spikes or microstructure-related disturbances.
The likelihood theory developed here may also be useful for modern machine-learning approaches to Voigt fitting. Recent work in MRI spectroscopy and quasar absorption-line analysis uses machine learning to accelerate Voigt-profile fitting by training algorithms to reproduce least-squares or traditional profile-fitting outputs. Our results complement this approach by characterizing the likelihood-based statistical target itself. In particular, the score and Fisher information can be used to assess local parameter identification, construct uncertainty measures, and provide one-step likelihood corrections to fast machine-learning estimates.
Overall, the paper shows that the Gauss-Cauchy convolution is not merely a useful heavy-tailed density, but a tractable framework for likelihood inference, robust signal extraction, and state-space filtering. Natural extensions include multivariate GCC models, richer latent-state dynamics, hybrid machine-learning and likelihood-based estimation, and applications in other settings where persistent Gaussian variation is observed through intermittent heavy-tailed measurement noise.
\bibliographystyle{apalike}
\bibliography{prh}