Extracted main text — title through conclusion, appendix excluded. This is what our citation measures are computed over, published so the extraction can be checked by eye.
48,472 characters · 7 sections · 27 citation commands
Dynamic Spatiotemporal ARCH Models: Small and Large Sample Results
\sloppy
\singlespacing
JEL-Classification: C11, C23, C58.\\ Keywords: Spatial ARCH, GMM, QMLE, volatility clustering, volatility, spatial dependence. \onehalfspacing
This paper investigates the small and large sample properties of three estimators for dynamic spatiotemporal ARCH models suggested by Otto:2023. This model allows the log-volatility term to depend on (i) the spatial lag of the log-squared outcome variable, (ii) the time-lag of the log-squared outcome variable, (iii) the spatiotemporal lag of the log-squared outcome variable, (iv) exogenous variables, and (v) the unobserved heterogeneity across regions and time, i.e., the regional and time fixed effects. The estimation equation of the model is obtained through a log-squared transformation Robinson:2009, Taspinar:2021. Notably, while the estimation equation obtained is in the form of a standard spatial dynamic panel data model considered by Yu:2008 and Lee:2010, it incorporates two important new features.
Firstly, the outcome variable, the spatial, temporal, and spatiotemporal lags are formulated using the log-squared original outcome variable. Secondly, the disturbance term in the model is the log-square of the original disturbance term due to the log-squared transformation. Following Lee:2007 and Lee:2014, Otto:2023 consider a generalized method of moment (GMM) estimator based on a set of linear and quadratic moment functions. In this paper, we also consider two quasi-maximum likelihood (QML) estimators considered in Yu:2008 and Lee:2010 for the estimation of the model. The first estimator is based on a transformation approach requiring the estimation of regional fixed effects along with the other model parameters. The second estimator is based on a direct approach and necessitates the estimation of both regional and time fixed effects. We first compare the theoretical properties of these estimators and subsequently investigate their small-sample properties through Monte Carlo simulations.
The rest of the paper proceeds in the following way. In Section (ref), we describe the dynamic spatiotemporal ARCH model and show how it differs from a standard spatial dynamic panel data model. In Section (ref), we summarize the QML and the GMM estimation approaches and provide their large sample results. In Section (ref), we investigate the finite sample properties of both estimation approaches through Monte Carlo simulations. In Section (ref), we conclude and provide some directions for future studies.
We consider the random process $\{Y_t(\mathbf{s}): \mathbf{s}\in \mathbf{D}_1\subseteq\mathbb{R}^d, d>1, t\in D_2\subseteq\mathbb{R}\}$, where $\mathbf{s}\in \mathbf{D}_1$ denotes the spatial location, and $t\in D_2$ is the time point. The structure of the spatial domain $\mathbf{D}_1$ and the time domain $D_2$ depends on the nature of spatial data, and we will assume that $\mathbf{D}_1=\{\mathbf{s}_1,\hdots,\mathbf{s}_n\}$ and $D_2=\{1,2,\hdots,T\}$. Then, the dynamic spatiotemporal ARCH model suggested by Otto:2023 can be expressed as
where $h_{t}(\mathbf{s}_i)$ is considered as the volatility term in location $\mathbf{s}_i$ at time $t$, and $\operatorname{\varepsilon}_{it}(\mathbf{s}_i)$ are independent and identically distributed random variables that have mean zero and unit variance. The log-volatility terms follow the process specified in (ref), where $\{m_{ij}\}$, for $i,j=1,\hdots,n$, are the non-stochastic spatial weights, $\mathbf{x}_{it}$ is a $k\times1$ vector of exogenous variables with the associated parameter vector $\boldsymbol{\beta}_0$, and $\boldsymbol{\mu}_0=(\mu_{0}(\mathbf{s}_1),\hdots,\mu_{0}(\mathbf{s}_n))^{'}$ and $\boldsymbol{\alpha}_0=(\alpha_{10},\hdots,\alpha_{T0})^{'}$ are spatial and time fixed effects. In the log-volatility equation, the spatial, temporal and spatiotemporal effects of the log-squared outcome variable on the log-volatility are measured by the unknown scalar parameters $\gamma_0$, $\rho_{0}$, and $\delta_{0}$, respectively. We assume that both $\boldsymbol{\mu}_0$ and $\boldsymbol{\alpha}_0$ can be correlated with the exogenous variables in an arbitrary manner and the initial value vector $\mathbf{Y}_0=(Y_{0}(\mathbf{s}_1),\hdots,Y_{0}(\mathbf{s}_n))^{'}$ is observable.
Define $Y^*_{t}(\mathbf{s}_i)=\log Y^2_{t}(\mathbf{s}_i)$, $h^{*}_{t}(\mathbf{s}_i)=\log h_{t}(\mathbf{s}_i)$ and $\operatorname{\varepsilon}^*_{t}(\mathbf{s}_i)=\log\operatorname{\varepsilon}^2_{t}(\mathbf{s}_i)$. Then, we apply the log-squared transformation to (ref) and obtain
In vector form, we can express (ref) and (ref) in the following way:
where $\mathbf{M}=(m_{ij})$ is the $n\times n$ matrix of the spatial weights, $\mathbf{Y}^{*}_t=(Y^{*}_{t}(\mathbf{s}_1),\hdots,Y^{*}_{t}(\mathbf{s}_n))^{'}$, $\mathbf{h}^{*}_t=(h^{*}_{t}(\mathbf{s}_1),\hdots,h^{*}_{t}(\mathbf{s}_n))^{'}$, $\boldsymbol{\operatorname{\varepsilon}}^{*}_t=(\operatorname{\varepsilon}^{*}_{t}(\mathbf{s}_1),\hdots,\operatorname{\varepsilon}^{*}_{t}(\mathbf{s}_n))^{'}$, $\mathbf{X}_t=(\mathbf{x}_{t}(\mathbf{s}_1),\hdots,\mathbf{x}_{t}(\mathbf{s}_n))^{'}$, and $\mathbf{1}_n$ is the $n\times1$ vector of ones. Then, we obtain the following estimation equation by substituting (ref) into (ref):
The estimation equation in (ref) is in the form of a spatial dynamic panel data model considered by Yu:2008 and Lee:2010. However, it differs from a standard spatial dynamic panel data model in two important ways, which have implications for the chosen estimation approach. Firstly, the outcome variable, the spatial lag term ($\mathbf{M}\mathbf{Y}^{*}_t$), and the spatiotemporal lag term ($\mathbf{M}\mathbf{Y}^{*}_{t-1}$) are formulated in terms of the log-squared outcome variable. Secondly, the elements of $\boldsymbol{\operatorname{\varepsilon}}^{*}_t$ are the log-squared original disturbance terms, i.e., $\operatorname{\varepsilon}^*_{t}(\mathbf{s}_i)=\log\operatorname{\varepsilon}^2_{t}(\mathbf{s}_i)$. If we assume that $\operatorname{\varepsilon}_{t}(\mathbf{s}_i)\sim N(0,1)$, then $\operatorname{\varepsilon}^{*}_{t}(\mathbf{s}_i)\sim\log\chi^2_1$, which is the log-chi squared distribution with one degree of freedom with the following density:
The first two moments of $\log\chi^2_1$ are $\mathbb{E}(\operatorname{\varepsilon}^{*}_{t}(\mathbf{s}_i)) = -\- \log(2) \approx-1.2704$, where $\gamma$ is Euler's constant, and $\text{Var}(\operatorname{\varepsilon}^{*}_{t}(\mathbf{s}_i))=\pi^2/2\approx4.9348$ Peter:2012. In Figure (ref), we compare this density with that of $N(-\- \log(2),\,\pi^2/2)$. As seen from the figure, the log-chi squared distribution exhibits significant skewness with a long left tail.
Although the elements of $\boldsymbol{\operatorname{\varepsilon}}^{*}_t$ in (ref) are i.i.d across $\mathbf{s}_i$ and $t$, we may have $\mathbb{E}\left(\boldsymbol{\operatorname{\varepsilon}}^{*}_t\right)\ne\mathbf{0}$ because of the log-squared transformation. Let $\mathbb{E}\left(\boldsymbol{\operatorname{\varepsilon}}^{*}_t\right)=\mu_{\operatorname{\varepsilon}}\mathbf{1}_n$, where $\mu_{\operatorname{\varepsilon}}$ is a scalar unknown parameter, and define $\mathbf{U}_t=(u_{t}(\mathbf{s}_1),\hdots,u_{t}(\mathbf{s}_n))^{'}=\boldsymbol{\operatorname{\varepsilon}}^{*}_t-\mathbb{E}\left(\boldsymbol{\operatorname{\varepsilon}}^{*}_t\right)=\boldsymbol{\operatorname{\varepsilon}}^{*}_t-\mu_{\operatorname{\varepsilon}}\mathbf{1}_n$. Then, we can express the estimation equation as
Following the spatial econometric literature, we consider two estimation approaches for (ref): (i) the QML method, and (ii) the GMM method. In both approaches, we assume that $\{u_{t}(\mathbf{s}_i)\}$ are i.i.d across $t$ and $\mathbf{s}_i$ with mean zero and variance $\sigma^2_0$.
Following Lee:2010, we consider two QML methods: (i) the transformation approach and (ii) the direct approach. In the case of the transformation approach, we will eliminate the time fixed effects from the model through a suitable transformation and then estimate the regional fixed effects along with the other parameters. In the direct approach, we will estimate both $\boldsymbol{\alpha}_0$ and $\mathbf{c}_{n0}=\boldsymbol{\mu}_0+\mu_{\operatorname{\varepsilon}}\mathbf{1}_n$ simultaneously. The transformation approach is only applicable if $\mathbf{M}$ is row-normalized, i.e., if $\mathbf{M}\mathbf{1}_n=\mathbf{1}_n$, while the direct approach does not require the row-normalization.
We start with the transformation approach and assume that $\mathbf{M}\mathbf{1}_n=\mathbf{1}_n$. Let $\mathbf{J}_n=\mathbf{I}_n-\frac{1}{n}\mathbf{1}_n\mathbf{1}^{'}_n$, where $\mathbf{I}_n$ is the $n\times n$ the identity matrix, and $(\mathbf{F}_{n,n-1},\mathbf{1}_n/\sqrt{n})$ be the orthonormal matrix of eigenvectors of $\mathbf{J}_n$, where $\mathbf{F}_{n,n-1}$ is the $n\times (n-1)$ matrix of eigenvectors corresponding to the eigenvalues of ones.\footnote{Some properties of $\mathbf{F}_{n,n-1}$ are $\mathbf{J}_n\mathbf{F}_{n,n-1}=\mathbf{F}_{n,n-1}$, $\mathbf{F}^{'}_{n,n-1}\mathbf{F}_{n,n-1}=\mathbf{I}_{n-1}$, $\mathbf{F}^{'}_{n,n-1}\mathbf{1}_{n}=\mathbf{0}$, $\mathbf{F}_{n,n-1}\mathbf{F}^{'}_{n,n-1}+\frac{1}{n}\mathbf{1}_n\mathbf{1}^{'}_n=\mathbf{I}_n$ and $\mathbf{F}_{n,n-1}\mathbf{F}^{'}_{n,n-1}=\mathbf{I}_n$. } Then, since $\mathbf{F}^{'}_{n,n-1}\mathbf{1}_n=\mathbf{0}$ holds, multiplying both sides of (ref) with $\mathbf{F}^{'}_{n,n-1}$ yields
where $\mathbf{Y}^{**}_t=\mathbf{F}^{'}_{n,n-1}\mathbf{Y}^{*}_t$, $\mathbf{M}^{*}=\mathbf{F}^{'}_{n,n-1}\mathbf{M}\mathbf{F}_{n,n-1}$, $\mathbf{X}^{*}_t=\mathbf{F}^{'}_{n,n-1}\mathbf{X}_t$, $\boldsymbol{\mu}^{*}_0=\mathbf{F}^{'}_{n,n-1}\boldsymbol{\mu}_0$ and $\mathbf{U}^{*}_t=\mathbf{F}^{'}_{n,n-1}\mathbf{U}_t$. Let $\boldsymbol{\theta}_0=(\gamma_0,\delta_0,\boldsymbol{\beta}^{'}_0,\rho_0,\sigma^2_0)^{'}$ denote the true parameter vector and $\boldsymbol{\theta}=(\gamma,\delta,\boldsymbol{\beta}^{'},\rho,\sigma^2)^{'}$ denote any other arbitrary value. If we assume that $\mathbf{U}_t\sim N(\mathbf{0},\sigma^2_0\mathbf{I}_n)$, then we have $\mathbf{U}^{*}_t\sim N(\mathbf{0},\sigma^2_0\mathbf{I}_{n-1})$. Thus, the log-likelihood function of (ref) can be expressed as
where $\mathbf{U}^{*}_t(\boldsymbol{\theta})=(\mathbf{I}_{n-1}-\rho\mathbf{M}^{*})\mathbf{Y}^{**}_{t}-\mathbf{Z}^{**}_t\boldsymbol{\eta}-\boldsymbol{\mu}^{*}$ with $\mathbf{Z}^{**}_t=(\mathbf{Y}^{**}_{t-1},\mathbf{M}^{*}\mathbf{Y}^{**}_{t-1},\mathbf{X}^{*}_t)$ and $\boldsymbol{\eta}=(\gamma,\delta,\boldsymbol{\beta}^{'})^{'}$. Using the properties of $\mathbf{F}_{n,n-1}$, we can show that (i) $\ln|\mathbf{I}_{n-1}-\rho\mathbf{M}^{*}|=1/(1-\rho)|\mathbf{I}_n-\rho\mathbf{M}|$, and (ii) $\sum_{t=1}^T\mathbf{U}^{*'}_t(\boldsymbol{\theta})\mathbf{U}^{*}_t(\boldsymbol{\theta})=\sum_{t=1}^T\mathbf{U}^{'}_t(\boldsymbol{\theta})\mathbf{J}_n\mathbf{U}_t(\boldsymbol{\theta})$, where $\mathbf{U}_t(\boldsymbol{\theta})=(\mathbf{I}_{n-1}-\rho\mathbf{M})\mathbf{Y}^{*}_{t}-\mathbf{Z}^{*}_t\boldsymbol{\eta}-\boldsymbol{\mu}$ with $\mathbf{Z}^{*}_t=(\mathbf{Y}^{*}_{t-1},\mathbf{M}\mathbf{Y}^{*}_{t-1},\mathbf{X}_t)$. Thus, we can express (ref) in terms of the original variables as
For an $n\times1$ vector $\mathbf{V}_t$, we define $\tilde{\mathbf{V}}_t=\mathbf{V}_t-\frac{1}{T}\sum_{t=1}^T\mathbf{V}_t$. Then, concentrating out $\boldsymbol{\mu}$ from (ref) yields
where $\tilde{\mathbf{U}}_t(\boldsymbol{\theta})=(\mathbf{I}_{n-1}-\rho\mathbf{M})\tilde{\mathbf{Y}}^{*}_{t}-\tilde{\mathbf{Z}}^{*}_t\boldsymbol{\eta}-\boldsymbol{\mu}$ with $\tilde{\mathbf{Z}}^{*}_t=(\tilde{\mathbf{Y}}^{*}_{t-1},\mathbf{M}\tilde{\mathbf{Y}}^{*}_{t-1},\tilde{\mathbf{X}}_t)$. Then, the QMLE $\widehat{\text{$\boldsymbol{\theta}$}}_{nT}$ of $\boldsymbol{\theta}_0$ is defined by $\widehat{\text{$\boldsymbol{\theta}$}}_{nT}=\operatorname{argmax}_{\boldsymbol{\theta}}\ln L(\boldsymbol{\theta})$. To investigate the asymptotic properties, following Lee:2010, we can show that the score functions can be decomposed as $\frac{1}{\sqrt{(n-1)T}}\frac{\partial\ln L(\boldsymbol{\theta}_0)}{\partial\boldsymbol{\theta}}=\frac{1}{\sqrt{(n-1)T}}\frac{\partial\ln L^u(\boldsymbol{\theta}_0)}{\partial\boldsymbol{\theta}}-\Delta_{nT}$, where the first component is uncorrelated with $\mathbf{U}_t$ while the second component is correlated with $\mathbf{U}_t$ for $t\leq T-1$. Moreover, the second component is the source of the asymptotic bias with the order $\Delta_{nT}=\sqrt{\frac{n-1}{T}}\mathbf{a}+O(\sqrt{(n-1)/T^3})+O_p(1/T^3)$, where $\mathbf{a}=O(1)$. Then, following Lee:2010, it can be shown that
where $\boldsymbol{\Sigma}=\lim_{T\to\infty}\mathbb{E}\left(-\frac{1}{(n-1)T}\frac{\partial^2\ln L(\boldsymbol{\theta}_0)}{\partial\boldsymbol{\theta}\partial\boldsymbol{\theta}^{'}}\right)$, $\boldsymbol{\Omega}=\lim_{T\to\infty}\text{Var}\left(\frac{1}{\sqrt{(n-1)T}}\frac{\partial\ln L^u(\boldsymbol{\theta}_0)}{\partial\boldsymbol{\theta}}\right)-\boldsymbol{\Sigma}$, and $\mathbf{b}=\boldsymbol{\Sigma}^{-1}\mathbf{a}=O(1)$ is the asymptotic bias term.\footnote{The explicit forms of $\mathbf{b}$, $\boldsymbol{\Sigma}$ and $\boldsymbol{\Omega}$ can be readily determined using the results presented in Lee:2010. } According to the asymptotic result in (ref), there are three cases depending on the relative growth rates of $n$ and $T$. Firstly, if $\frac{n}{T}\to\infty$, i.e., if $n$ grows faster than $T$, then $T(\widehat{\text{$\boldsymbol{\theta}$}}_{nT}-\boldsymbol{\theta}_0)+\mathbf{b}\xrightarrow{p}\mathbf{0}$. That is, $\widehat{\text{$\boldsymbol{\theta}$}}_{nT}$ is consistent with rate $T$ and has a degenerate limiting distribution. Secondly, if $\frac{n}{T}\to0$, then we will have $\sqrt{T(n-1)}(\widehat{\text{$\boldsymbol{\theta}$}}_{nT}-\boldsymbol{\theta}_0)\xrightarrow{d}N\left(\mathbf{0},\,\boldsymbol{\Sigma}^{-1}(\boldsymbol{\Sigma}+\boldsymbol{\Omega})\boldsymbol{\Sigma}^{-1}\right)$. Lastly, when $\frac{n}{T}\to c\in \mathbb{R}_+$, i.e., when $T$ is asymptotically proportional to $n$, we have $\sqrt{T(n-1)}(\widehat{\text{$\boldsymbol{\theta}$}}_{nT}-\boldsymbol{\theta}_0)+\sqrt{c}\mathbf{b}\xrightarrow{d}N\left(\mathbf{0},\,\boldsymbol{\Sigma}^{-1}(\boldsymbol{\Sigma}+\boldsymbol{\Omega})\boldsymbol{\Sigma}^{-1}\right)$. In this case, Lee:2010 suggest a bias corrected estimator defined by $\widehat{\text{$\boldsymbol{\theta}$}}^1_{nT}=\widehat{\text{$\boldsymbol{\theta}$}}_{nT}-\widehat{\text{$\mathbf{b}$}}$, where $\widehat{\text{$\mathbf{b}$}}$ is a plug-in estimator of $\mathbf{b}$. Then, it follows that $\sqrt{T(n-1)}(\widehat{\text{$\boldsymbol{\theta}$}}^1_{nT}-\boldsymbol{\theta}_0)\xrightarrow{d}N\left(\mathbf{0},\,\boldsymbol{\Sigma}^{-1}(\boldsymbol{\Sigma}+\boldsymbol{\Omega})\boldsymbol{\Sigma}^{-1}\right)$ when $\frac{n}{T^3}\to0$, i.e., when $T^3$ grows faster than $n$.
Next, we introduce the direct approach for the estimation of (ref). Under the assumption that $\mathbf{U}_t\sim N(\mathbf{0},\sigma^2_0\mathbf{I}_n)$, we can derive the log-likelihood function of (ref) as
where $\mathbf{U}_t(\boldsymbol{\theta},\mathbf{c}_n,\boldsymbol{\alpha})=(\mathbf{I}_{n}-\rho\mathbf{M})\mathbf{Y}^{*}_{t}-\mathbf{Z}^{*}_t\boldsymbol{\eta}-\mathbf{c}_n+\alpha_t\mathbf{1}_n$. Concentrating both $\mathbf{c}_n$ and $\boldsymbol{\alpha}$ from (ref) yields
The QMLE $\widehat{\text{$\boldsymbol{\theta}$}}^d_{nT}$ of $\boldsymbol{\theta}_0$ is defined by $\widehat{\text{$\boldsymbol{\theta}$}}^d_{nT}=\operatorname{argmax}_{\boldsymbol{\theta}}\ln L^d(\boldsymbol{\theta})$. As in the case of the transformation approach, the score functions can be decomposed as $\frac{1}{\sqrt{nT}}\frac{\partial \ln L^d(\boldsymbol{\theta})}{\partial\boldsymbol{\theta}}=\frac{1}{\sqrt{nT}}\frac{\partial \ln L^{du}(\boldsymbol{\theta})}{\partial\boldsymbol{\theta}}-\Delta_{1,nT}-\Delta_{2,nT}$, where the first component is uncorrelated with $\mathbf{U}_t$ and the other components are the sources of asymptotic bias with the following orders: $\Delta_{1,nT}=\sqrt{\frac{n}{T}}\mathbf{a}_1+O(\sqrt{n/T^3})+O_p(1/\sqrt{T})$ with $\mathbf{a}_1=O(1)$, and $\Delta_{2,nT}=\sqrt{\frac{n}{T}}\mathbf{a}_2$ with $\mathbf{a}_2=O(1)$. Then, it follows from Lee:2010 that
where $\boldsymbol{\Sigma}=\lim_{T\to\infty}\mathbb{E}\left(-\frac{1}{nT}\frac{\partial^2\ln L^d(\boldsymbol{\theta}_0)}{\partial\boldsymbol{\theta}\partial\boldsymbol{\theta}^{'}}\right)$, $\boldsymbol{\Omega}=\lim_{T\to\infty}\text{Var}\left(\frac{1}{\sqrt{nT}}\frac{\partial\ln L^{du}(\boldsymbol{\theta}_0)}{\partial\boldsymbol{\theta}}\right)-\boldsymbol{\Sigma}$, and $\mathbf{b}_1=\boldsymbol{\Sigma}^{-1}\mathbf{a}_1=O(1)$ and $\mathbf{b}_2=\boldsymbol{\Sigma}^{-1}\mathbf{a}_2=O(1)$ are the asymptotic bias terms. As in the transformation case, there are three cases depending on the relative growth rates of $n$ and $T$. The first two cases are the degenerate limiting distribution cases, occurring when either $\frac{n}{T}\to0$ or $\frac{n}{T}\to\infty$ holds. If $\frac{n}{T}\to0$ holds, then we have $n(\widehat{\text{$\boldsymbol{\theta}$}}^d_{nT}-\boldsymbol{\theta}_0)+\mathbf{b}_2\xrightarrow{p}\mathbf{0}$, and when $\frac{n}{T}\to\infty$ holds, we have $T(\widehat{\text{$\boldsymbol{\theta}$}}^d_{nT}-\boldsymbol{\theta}_0)+\mathbf{b}_1\xrightarrow{p}\mathbf{0}$. Finally, when $\frac{n}{T}\to c\in \mathbb{R}_+$ holds, we have $\sqrt{nT)}(\widehat{\text{$\boldsymbol{\theta}$}}^d_{nT}-\boldsymbol{\theta}_0)+\sqrt{c}\mathbf{b}_1+\sqrt{1/c}\mathbf{b}_2\xrightarrow{d}N\left(\mathbf{0},\,\boldsymbol{\Sigma}^{-1}(\boldsymbol{\Sigma}+\boldsymbol{\Omega})\boldsymbol{\Sigma}^{-1}\right)$, which indicates that we need a bias correction for $\widehat{\text{$\boldsymbol{\theta}$}}^d_{nT}$. Let $\boldsymbol{\theta}^{d1}_{nT}=\widehat{\text{$\boldsymbol{\theta}$}}^d_{nT}-\widehat{\text{$\mathbf{b}$}}_1/T-\widehat{\text{$\mathbf{b}$}}_2/T$ be the bias corrected estimator, where $\widehat{\text{$\mathbf{b}$}}_1$ and $\widehat{\text{$\mathbf{b}$}}_2$ are the plug-in estimators of $\mathbf{b}_1$ and $\mathbf{b}_2$, respectively. Then, it can be shown that if $\frac{n}{T^3}\to0$ and $\frac{T}{n^3}\to\infty$, then we have $\sqrt{nT)}(\widehat{\text{$\boldsymbol{\theta}$}}^{d1}_{nT}-\boldsymbol{\theta}_0)\xrightarrow{d}N\left(\mathbf{0},\,\boldsymbol{\Sigma}^{-1}(\boldsymbol{\Sigma}+\boldsymbol{\Omega})\boldsymbol{\Sigma}^{-1}\right)$.
In the GMM approach, we will use two different transformations to eliminate $\boldsymbol{\mu}$, $\boldsymbol{\alpha}$ and $\mu_{\operatorname{\varepsilon}}$ from the model. The first transformation is based on the decomposition of $\mathbf{J}_T=\left(\mathbf{I}_T-\frac{1}{T}\mathbf{1}_T\mathbf{1}^{'}_T\right)$, where $\mathbf{I}_T$ is the $T\times T$ identity matrix. Let $\left(\mathbf{F}_{T,T-1},\mathbf{1}_T/\sqrt{T}\right)$ be the orthonormal eigenvector matrix of $\mathbf{J}_T$, where $\mathbf{F}_{T,T-1}$ is the $T\times (T-1)$ sub-matrix containing eigenvectors corresponding to the eigenvalues of one. Let $\mathbf{D}=\left(\mathbf{d}_1,\hdots,\mathbf{d}_T\right)$ be an $n\times T$ matrix, where $\mathbf{d}_t$ is an $n\times1$ column vector for $t=1,\hdots,T$. Using $\mathbf{F}_{T,T-1}$, we can transform $\mathbf{D}$ into an $n\times(T-1)$ matrix in the following way: $\mathbf{D}^{*}=\left(\mathbf{d}^{*}_1,\hdots,\mathbf{d}^{*}_{T-1}\right)=\left(\mathbf{d}_1,\hdots,\mathbf{d}_T\right)\mathbf{F}_{T,T-1}$. Since $\left(\boldsymbol{\mu}_0,\hdots,\boldsymbol{\mu}_0\right)\mathbf{F}_{T,T-1}=\boldsymbol{\mu}_0\mathbf{1}^{'}_T\mathbf{F}_{T,T-1}=\mathbf{0}_{n\times(T-1)}$ and $\left(\mu_{\operatorname{\varepsilon}}\mathbf{1}_n,\hdots,\mu_{\operatorname{\varepsilon}}\mathbf{1}_n\right)\mathbf{F}_{T,T-1}=\mu_{\operatorname{\varepsilon}}\mathbf{1}_n\mathbf{1}^{'}_T\mathbf{F}_{T,T-1}=\mathbf{0}_{n\times(T-1)}$, this transformation can remove both $\boldsymbol{\mu}_0$ and $\mu_{\operatorname{\varepsilon}}$ from the model. If we apply $\mathbf{F}_{T,T-1}$ to (ref) in a similar manner, we will obtain
Note that (ref) still includes the transformed time fixed effects denoted by $\alpha^{*}_{t0}$ for $t=1,\hdots,T$. To remove these effects, we apply a second transformation by pre-multiplying (ref) with $\mathbf{J}_n$:
Following Lee:2007 and Lee:2014, Otto:2023 consider both linear and quadratic moment functions for the estimation of (ref). Let $N=n(T-1)$, $\boldsymbol{\theta}=(\rho,\boldsymbol{\eta}^{'})^{'}$ and $\mathbf{J}_N=\mathbf{I}_{T-1}\otimes\mathbf{J}_n$. Then, the vector of moment functions takes the following form:
where $\mathbf{U}_N(\boldsymbol{\theta})=(\mathbf{U}^{*'}_1(\boldsymbol{\theta}),\hdots,\mathbf{U}^{*'}_{T-1}(\boldsymbol{\theta}))^{'}$ with $\mathbf{U}^{*}_t(\boldsymbol{\theta})=(\mathbf{I}_n-\rho\mathbf{M})\mathbf{Y}^{**}_t-\mathbf{Z}^{**}_t\boldsymbol{\eta}-\alpha^{*}_{t}\mathbf{1}_n$. In (ref), the linear moment function takes the form of $\mathbf{Q}^{'}_N\mathbf{J}_N\mathbf{U}_N(\boldsymbol{\theta})$, where $\mathbf{Q}_N$ is the $n(T-1)\times k_q$ matrix of IVs, and the quadratic moment function takes the form of $\mathbf{U}^{'}_N(\boldsymbol{\theta})\mathbf{J}_N\mathbf{P}_{jN}\mathbf{J}_N\mathbf{U}_N(\boldsymbol{\theta})$, where $\mathbf{P}_{jN}=\mathbf{I}_{T-1}\otimes\mathbf{P}_j$ and $\mathbf{P}_j$ is an $n\times n$ matrix satisfying $\mathrm{tr}(\mathbf{J}_n\mathbf{P}_j\mathbf{J}_n)=0$ for $j=1,2,\hdots,m$.
Let $\boldsymbol{\Omega}_N=\frac{1}{N}\mathbb{E}\left(\mathbf{g}_N(\boldsymbol{\theta}_0)\mathbf{g}^{'}_N(\boldsymbol{\theta}_0)\right)$. Then, the optimal GMME is defined by
where $\widehat{\text{$\boldsymbol{\Omega}$}}_N$ be a consistent estimator of $\boldsymbol{\Omega}_N$, i.e., $\widehat{\text{$\boldsymbol{\Omega}$}}_N-\boldsymbol{\Omega}_N=o_p(1)$. Under some assumptions, Otto:2023 show that $\frac{1}{N}\frac{\partial \mathbf{g}_N(\boldsymbol{\theta}_0)}{\partial\boldsymbol{\theta}^{'}}=\mathbf{D}_{1N}+\mathbf{D}_{2N}+O_p(N^{-1/2})$, where $\mathbf{D}_{1N}=O(1)$ and $\mathbf{D}_{2N}=O(T^{-1})$.\footnote{See Otto:2023 for the explicit forms of $\mathbf{D}_{1N}$ and $\mathbf{D}_{2N}$.} Then, when $T$ is finite and $n\to\infty$, it can be shown that
On the other hand, when $T\to\infty$, we have $\mathbf{D}_{2N}=O(T^{-1})=o(1)$. Thus, when $T\to\infty$ and $n\to\infty$, Otto:2023 show that
Note that the asymptotic results in (ref) and (ref) hold for the arbitrary IV matrix $\mathbf{Q}_N$ and the quadratic moment matrices $\mathbf{P}_{jN}$ for $j=1,2,\hdots,m$. The asymptotic efficiency of $\widehat{\text{$\boldsymbol{\theta}$}}_N$ should be considered in choosing the IV and quadratic moment matrices. Following Lee:2014, Otto:2023 provide a vector of moment functions that can lead to the most efficient GMME when both $n$ and $T$ are large.
In this section, we investigate the finite sample performance of all three considered estimators when the sample size is small. In particular, we are looking at the case of short time series (i.e., $T$ is small), which is often encountered in practice for geo-referenced data. More precisely, we simulate $Y_{t}(\boldsymbol{s}_i) = h_{t}(\boldsymbol{s}_i)^{1/2}\varepsilon_{t}(\boldsymbol{s}_i)$ with
We consider three different data-generating processes with the following parameter settings:
In $M_1$, the temporal effect is relatively dominant; in $M_2$, all effects are relatively weak; and finally, in $M_3$, the spatial effect is relatively strong. Both $M_1$ and $M_2$ include only the spatial fixed effects, while $M_3$ includes both the spatial and time fixed effects. The spatial/time fixed effects and the exogenous regressors are independently simulated from $N(0,1)$. We generate the error terms $\varepsilon_{t}(\boldsymbol{s}_i)$'s independently from $N(0,1)$, and consider a queen contiguity row-normalized weights matrix. The sample size gradually increases for both $n$ and $T$ with $n \in \{25, 49, 81\}$ on a regular two-dimensional lattice (with side length 5, 7, and 9) and $T \in \{ 5, 10, 20\}$. The number of repetitions is set to $1000$ in all cases. Thus, we considered 27 different model specifications in total and three different estimators, which were always applied to the same simulated values.
For each estimator, we report bias (BIAS), the root mean square error (RMSE) and the mean absolute error (MAE). The results are reported in Tables (ref) and (ref).\footnote{The MAE results are provided in the accompanying web appendix. These results suggest the same conclusions as those obtained from the RMSE results presented in Table (ref). Additionally, in the accompanying web appendix, we provide graphical representations of the results provided in Tables (ref) and (ref).}
In Table (ref), we observe that all estimators report relatively larger bias when $(n,T)=(25, 5)$, especially in the case of $\rho$, $\gamma$ and $\delta$. When either $n$ or $T$ increases, all estimators impose a smaller bias in all cases. In the case of $\rho$, the QMLE based on the direct approach performs relatively better than other estimators. In the case of $\gamma$, the GMME performs better than the QML estimators. It is clear that the QML methods require a large $T$ to estimate this parameter accurately. In the case of $\delta$, the GMME seems to impose a relatively smaller bias than the QML estimators. For $\beta_0$ and $\beta_1$, the QMLE based on the transformation approach shows a significant bias when $T=5$.
In Table (ref), we observe that the RMSE decreases when either $n$ or $T$ increases for all estimators. In the case of $\rho$, the QMLE based on the direct approach is relatively more efficient than the other estimators. For instance, in $M_2$ with $(n,T)=(25,20)$, the RMSEs are $0.185$, $0.148$ and $0.063$ for the GMME, the QMLE based on transformation and the QMLE based on the direct approach, respectively. In the case of $\gamma$, the GMME is relatively more efficient than the QML estimators in all cases. In the case of $\delta$, the direct approach performs slightly better than the other approaches. All estimators perform similarly in the case of $\beta_0$ and $\beta_1$.
In this paper, we considered the small and large sample properties of (i) the QMLE based on the transformation approach, (ii) the QMLE based on the direct approach, and (iii) the GMME for the estimation of a dynamic spatiotemporal ARCH model. Theoretically, the estimators based on the QML method require a large $T$ and a bias correction approach. The bias-corrected QMLE based on the transformation approach requires that $\frac{n}{T^3}\to0$, while the bias-corrected QMLE based on the direct approach requires both $\frac{n}{T^3}\to0$ and $\frac{T}{n^3}\to\infty$. The QMLE based on the transformation approach is only applicable for models with row-normalized weights matrices. The GMME does not require bias correction and is valid under both finite and large $T$ cases. In a Monte Carlo simulation study, we compare the finite sample properties of these estimators. Our results indicate that the QMLE based on the direct approach performs relatively better for the estimation of the spatial effect ($\rho$), while the GMME performs relatively better for the estimation of the temporal effect ($\gamma$). Overall, when both $T$ and $n$ are large enough, all estimators perform similarly and satisfactorily. In future studies, the performance of these estimators can be investigated for a dynamic spatiotemporal ARCH model that has a multiplicative structure for the regional and time-fixed effects.