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.
84,309 characters · 17 sections · 84 citation commands
Functional Regression with Nonstationarity and Error Contamination: Application to the Economic Impact of Climate Change
\def\spacingset#1{ {#1}} \spacingset{1}
{\it Keywords:} Functional linear model, cointegration, measurement errors, climate change.
\onehalfspacing
In data-rich environments, practitioners often need to deal with non-traditional observations, such as curves, probability density functions, or images. Accordingly, recent literature on functional data analysis, which provides statistical methods for handling such complex data, has gained popularity. For a comprehensive and broad review of this topic, readers are referred to Ramsay2005 and HK2012. Practitioners in various fields have benefited from advances in this area. In particular, functional linear regression models have become a central tool for those interested in analyzing the relationships between two or more such variables. Some early contributions to this topic include Yao2005, Hall2007, Park2012397, Florence2015, Benatia2017 and imaizumi2018, and, more recently, Chen_et_al_2020, Babii2022 and seong2021functional study the issue of endogeneity. A common feature of all these papers is that they all consider functional regression models with iid or stationary sequence of random functions.
Only recently has the literature begun to consider nonstationary dependent observations, even if many economic and statistical functional time series tend to be nonstationary, as noted in recent papers (e.g., Chang2016152, BSS2017,Franchi2017b, LRS, NSS, NSS2, seo2020functional, seoshang22). As a result, statistical methods developed for such time series are currently limited to analyzing their essential properties, such as cointegration, stochastic trends, and the dominant subspace. Despite its empirical relevance, articles developing inferential methods for functional (auto-)regression models involving nonstationary functional time series are scarce; to the best of the authors' knowledge, there are currently only a few preprints (e.g., changshocks,Hu2016). We fill this gap by developing novel statistical methods for functional regression models where both regressand and regressor exhibit unit-root-type nonstationary behavior allowing cointegration, a feature particularly important for economic and financial applications.
In addition to incorporating nonstationarity into functional regression models, we aim to enhance the real-world applicability of our methods by addressing a typical and practical aspect of functional data that has recently been discussed in the literature: incomplete and partially observed data (see, e.g., Chen_et_al_2020; seong2021functional). In the majority of real data analyses, (i) each functional observation, say $x_t(u)$ for $u\in [a_1,a_2]$, is not directly observed, and often constructed from its partial and discrete realizations $(x_t(u_1),\ldots,x_t(u_n))'$ for $u_1,\ldots,u_n \in [a_1,a_2]$ and (ii) the number of discrete observations $n$ is often not large enough. In fact, (i) and (ii) are pointed out by seong2021functional in the context of functional linear models, and they argued that “endogeneity” caused by measurement errors need to be properly addressed for estimation and inference (see also Chen_et_al_2020). As a specific example, consider a case where functional observations are probability density-valued (as in Section (ref) to appear). This case has gained significant interest in the literature; for the stationary case, see e.g., kneip2001inference and Park2012397, while for the nonstationary case refer to Chang2016152 and seo2019cointegrated. In this scenario, the true probability density is not observable, and thus, it needs to be replaced by a proper nonparametric estimate. This naturally introduces small or large measurement errors in practice. In this paper, we explicitly consider cases where the variables of interest, which are nonstationary, are also error-contaminated, and then develop statistical methods that are robust to error contamination. This not only distinguishes the present paper significantly from existing works (cf., e.g., Benatia2017, Park2012397, Chen_et_al_2020, Babii2022,seong2021functional) but also makes our proposed methods more appealing to applied researchers. We also believe that our methodology can be applied to various economic and financial time series.
More technically and specifically, we assume that the variables of interest are cointegrated functional time series, following the framework of Chang2016152 and BSS2017. This assumption has been widely used in the recent literature on nonstationary functional time series, especially in economic applications (see, e.g., NSS, seo2020functional). We then assume that these variables can only be observed with additive measurement errors. As noted by seong2021functional, the problem of neglected error contamination generally results in inconsistency of standard estimators constructed from the sample covariance operator $\widehat{C}_0$ of the regressor (this is also true in our model, as will be discussed in Section (ref) in greater detail). This inconsistency arises primarily because $\widehat{C}_0$ is inherently contaminated by measurement errors and, consequently, becomes a distorted estimator of its population counterpart. To address this issue of error contamination, we consider autocovariance-based inference, avoiding the direct use of the covariance operator of an error-contaminated variable for statistical inference, as in some recent articles on functional regression models (see, e.g., Chen_et_al_2020). More specifically, we construct our proposed estimator based on the lag-$\kappa$ sample autocovariance $\widehat{C}_{\kappa}$ for some positive $\kappa$. This approach is grounded in the observation that, as long as the measurement errors are not strongly correlated and satisfy certain mild regularity conditions (to be detailed), (i) the sample autocovariance $\widehat{C}_{\kappa}$ will be less affected by measurement errors, and (ii) the assumption that measurement errors are not strongly correlated does not seem overly restrictive, given that such errors in functional data analysis commonly arise from constructing each individual functional observation based on its discrete realizations; as will be detailed in Section (ref), our asymptotic analysis requires a much weaker condition on the serial correlation of the measurement errors rather than complete serial uncorrelatedness. It should be noted that, in this paper, we also consider the case $\kappa=0$, which yields the standard covariance-based estimator in the functional linear model, and study its detailed asymptotic properties, a contribution that is, to the best of the authors' knowledge, also novel.
We develop autocovariance-based inferential methods that are robust to the potential presence of measurement errors. This includes a novel dimension-reduction method, our proposed estimator of the slope parameter in the functional regression model, and their asymptotic properties. The proposed estimator bears some resemblance to the conventional two-step estimator of engle1987 in that both use residuals computed from the estimated relationship between the nonstationary components of the model; however, beyond this superficial similarity, the two approaches differ substantially in structure and purpose. We also provide numerical studies with real-world data and simulation experiments to examine the performance of our proposed estimator. As an application, we illustrate the empirical relevance of our proposed methodology by applying it to an empirical model for studying the global economic impact of climate change. Specifically, we show that the proposed framework effectively estimates the distributional relationship between land temperature anomalies--often considered a measure of climate change--and regional economic growth rates in the possible presence of measurement errors, thereby offering a robust basis for assessing heterogeneous climate–economy relationships across the globe.
The rest of the paper is organized as follows. Section (ref) reviews essential preliminaries on nonstationary cointegrated functional time series. Section (ref) details the regression model, and Section (ref) develops inferential methods for it. In Section (ref), we apply the proposed method to examine the global economic impact of climate change. Section (ref) concludes.
We let \(\mathcal{H}\) be a real separable Hilbert space of functions on the interval $[a_1, a_2]$, and let $\langle \cdot, \cdot \rangle$ (resp.\ $\|\cdot\|$) denote the associated inner product (resp.\ norm). We let ${\mathcal H}_y$ denote another Hilbert space, which will be set to $\mathbb{R}$ (when the dependent variable $y_t$ is real-valued) or $\mathcal H$ (when $y_t$ is function-valued). Throughout, regardless of whether $\mathcal{H}_y = \mathbb{R}$ or $\mathcal{H}$, we adopt a slight abuse of notation by using $\langle \cdot, \cdot \rangle$ and $\|\cdot\|$ to denote the inner product and norm associated with $\mathcal{H}_y$, respectively. This notational simplification facilitates the exposition and poses minimal risk of confusion, as the meaning of each operation is readily inferred from the context. For the same reason, we use $I$ to denote the identity map on any Hilbert space under consideration. As a further simplification, we henceforth write $\int F$ to denote $\int_{0}^1 F(s)ds$ for any operator- or vector-valued function $F$ defined on $[0,1]$.
Section (ref) of the Supplementary Material reviews basic concepts of bounded linear operators and random elements associated with two (possibly different) Hilbert spaces. Accordingly, we let $\mathcal L_{\mathcal H}$ denote the space of bounded linear operators on $\mathcal H$ with the usual operator norm $\|\cdot\|_{\textrm{op}}$, and let $\otimes$ denote the tensor product associated with $\mathcal H$, $\mathcal H_y$, or both (see (ref)). Section (ref) also reviews $\mathcal H$-valued random elements $X$, their expectation (denoted $\mathbb{E}[X]$), covariance operator (denoted $C_X:=\mathbb{E}[(X-\mathbb{E}(X)) \otimes (X-\mathbb{E}(X))]$), and cross-covariance with an $\mathcal H_y$-valued random element $Y$ (denoted $C_{XY}:=\mathbb{E}[(X-\mathbb{E}(X)) \otimes (Y-\mathbb{E}(Y))]$). For $A \in \mathcal L_{\mathcal H}$, concepts such as the adjoint (denoted $A^\ast$), range (denoted $\operatorname{ran} A$), and kernel (denoted $\ker A$), as well as properties such as self-adjointness, compactness, Hilbert–Schmidtness, and nonnegativity are introduced in that section, and they will be useful subsequently.
We will consider sequences of random linear operators, constructed from random elements in $\mathcal H$ and $\mathcal H_y$ (for a more detailed discussion on general random linear operators, see skorohod2001). For any such operator-valued random sequence $\{A_j\}_{j\geq1}$, we write $A_j \to_{p} A$ to denote convergence in probability with respect to the operator norm (i.e., $\|A_j - A\|_{\textrm{op}} \to_p 0$). In the subsequent discussion, convergence in probability sometimes occurs for $\mathcal H$- or $\mathcal H_y$-valued elements (in the appropriate norm), but for convenience we use the same notation $\to_p$ to denote such convergence throughout, as distinguishing between the two would add notational complexity with little benefit. Moreover, as is common in the literature (see, e.g., seo2020functional), we write $A_j = A + O_p(a_T)$ (resp.\ $A_j = A + o_p(a_T)$) if $\|A_j - A\|_{\textrm{op}} = O_p(a_T)$ (resp.\ $\|A_j - A\|_{\textrm{op}} = o_p(a_T)$) for some sequence $a_T$. For any two operators $A$ and $B$, we write $A=_d B$ to denote equivalence in their finite-dimensional distributions as in seo2020functional, i.e., $A=_dB$ if, for any $n>0$, $\{v_j\}_{j=1}^n$ ($\subset \mathcal H \text{ or } \mathcal H_y$) and $\{w_j\}_{j=1}^n$ ($\subset \mathcal H \text{ or } \mathcal H_y$), the distribution of $(\langle A v_1,w_1 \rangle,\ldots,\langle A v_n,w_n \rangle)'$ equals that of $(\langle B v_1,w_1 \rangle,\ldots,\langle B v_n,w_n \rangle)'$.
We review cointegrated linear processes in $\mathcal H$, which have been used to model the persistent nonstationary behavior of many economic functional time series (see, e.g., Chang2016152, NSS, NSS2, seo2020functional). Suppose that $\Delta x_t = x_{t}-x_{t-1} = \sum_{j=0}^\infty \psi_j \varepsilon_{t-j}$ for some sequence of bounded linear operators $\{\psi_j\}_{j\geq 0}$ and an iid sequence $\{\varepsilon_t\}_{t \in \mathbb{Z}}$ satisfying $\mathbb{E}[\varepsilon_t] = 0$ and $\mathbb{E}[\|\varepsilon_t\|^4] < \infty$ and having a positive definite covariance $C_{\varepsilon}$. If $\sum_{j=0}^\infty j \|\psi_j\|_{\textrm{op}}<\infty$ holds, we know from the Phillips-Solo decomposition of Phillips1992 and its extension to a function space (see e.g., seo_2022), $x_t$ allows the following representation, ignoring the initial values that are negligible in our asymptotic analysis:
where $\psi(1) =\sum_{j=0}^\infty \psi_j$, $\eta_t = \sum_{j=0}^\infty \widetilde{\psi}_j \varepsilon_{t-j}$ and $\widetilde{\psi}_j = -\sum_{k=j+1}^\infty \psi_k$. Let ${\rm P}^S$ be the orthogonal projection onto $[\operatorname{ran} \psi(1)]^\perp$ and let ${\rm P}^N = I-{\rm P}^S$. Then $\langle x_t,v \rangle$ is stationary if and only if $v \in \mathcal H^S$ (see BSS2017). Thus the entire Hilbert space $\mathcal H$ can be orthogonally decomposed into $\mathcal H^N = \operatorname{ran} {\rm P}^N$ and $\mathcal H^S = \operatorname{ran} {\rm P}^S$. We call $\mathcal H^N$ (resp.\ $\mathcal H^S$) the nonstationary (resp.\ stationary) subspace induced by $\{x_t\}_{t\geq1}$.
Subsequently, we will consider cointegrated time series introduced in this section, but some additional restrictions will be imposed for our asymptotic analysis in Section (ref).
Let \(\{x_t\}_{t \geq 1}\) be a cointegrated \(\mathcal{H}\)-valued time series, as detailed in Section (ref), which induces a bipartite partition of \(\mathcal{H}\) into a nonstationary subspace \(\mathcal{H}^N\) and a stationary subspace \(\mathcal{H}^S\). We consider the following data-generating mechanism: for $f:\mathcal H \to {\mathcal H}_y$,
Note that the above model includes no deterministic terms. We first develop inferential methods for this case and then discuss extending our methods to a model with deterministic terms in Section (ref); as may be expected, this extension requires only modest and non-substantial modifications of the results developed for the case without deterministic terms.
Throughout this paper, we assume that $x_t$ cannot be directly observed but that only $\tilde{x}_t$, observed with measurement errors, is available. As highlighted by seong2021functional, this assumption is empirically relevant because functional observations used in practice are often incompletely observed, with only finitely many discrete realizations available to practitioners. Consequently, it is common to construct a functional observation $z_t$ in advance by smoothing its $n$ discrete data points $z_t(s_1), \ldots, z_t(s_n)$, with $s_j$ included in the entire interval $[a_1, a_2]$, before computing estimators or test statistics. While one may disregard measurement errors for simplicity if $n$ is large enough and the data points are densely observed over $[a_1, a_2]$, this is often not the case in practice (our empirical application in Section (ref) is an example). We develop theoretical results under the presence of measurement errors in functional variables, while also discussing how these results simplify in the absence of such errors. Accordingly, the subsequent theoretical developments remain applicable to the error-free case, which has been more commonly considered in the functional data analysis literature.
Particularly, the issue of measurement errors is prominent when considering probability density–valued functional observations, say $\{z_t\}_{t\geq1}$, or their relevant transformations $\{g(z_t)\}_{t\geq1}$ in practice. Since practitioners do not observe the true probability densities, they typically substitute them with appropriate nonparametric estimates in analysis, leading to inevitable estimation errors. As noted by seong2021functional, neglecting these estimation errors without proper treatment results in inconsistency of standard estimators used in functional linear models. In Section (ref), we consider a specific empirical example involving density-valued functional observations, providing a more detailed discussion based on the existing literature.
Specifically, we assume that $\tilde{x}_t$ is a measurement of $x_t$ with an additive error $e_t$, as follows:
where $e_t$ may generally be correlated with the variables $u_t^N$ and $u_t^S$ appearing (ref) and (ref), and it may also be serially correlated. A key assumption which we employ for our asymptotic analysis is that $e_{t-\kappa}$ (and also $e_{t+\kappa}$) for some finite $\kappa>0$ has asymptotically negligible sample (cross-)covariances with $u_t^N$, $u_t^S$, and $e_t$. Given that $e_{t-\kappa}$ or $e_{t+\kappa}$ is (primarily) smoothing error or random disturbance associated with a variable observed at a different time, this assumption is practically reasonable and likely satisfied even for small positive $\kappa$. Similarly, $y_t$ can suffer from error contamination, but as in conventional multivariate regression models, its measurement error is absorbed into $u_t$. This only changes the interpretation of $u_t$ in the subsequent analysis. In the presence of measurement errors $e_t$, we may rewrite (ref) as follows:
We introduce assumptions on the data-generating mechanism. Below, for notational convenience, we let $\widetilde{\mathcal H}= \mathcal H_y \times \mathcal H$ be the (Cartesian) product Hilbert space equipped with the inner product $\langle (h_1,h_2), (\ell_1,\ell_2) \rangle_{\widetilde{\mathcal H}} = \langle h_1,\ell_1 \rangle + \langle h_2,\ell_2 \rangle$ (note that we let $\langle \cdot,\cdot \rangle$ to denote the inner product on either $\mathcal H_y$ or $\mathcal H$ to simplify notation, so the former is the inner product on $\mathcal H_y$). Observing that $\mathcal H$ can be orthogonally decomposed by $\mathcal H^N$ and $\mathcal H^S$, we write $\mathcal H = \mathcal H^N \times \mathcal H^S$ and also write any $h\in \mathcal H$ as $({\rm P}^N h, {\rm P}^S h)$; of course, in this case, for any $(h_1,h_2) \in \mathcal H$ and $(\ell_1,\ell_2) \in \mathcal H$, the inner product on this product space can simply be represented by $\langle h_1,\ell_1 \rangle + \langle h_2,\ell_2 \rangle$ with the inner product associated with $\mathcal H$. We employ the following assumptions throughout: in the assumption below, we consider $\mathcal H$-valued process $\mathcal E_t^x = (u_t^N, u_t^S)$ and $\widetilde{\mathcal H}$-valued process $\mathcal E_t = (u_t, \mathcal E_t^x)$.
We also require assumptions on measurement errors. As detailed in Section (ref), our estimator relies on the lag-$\kappa$ autocovariance operator, with $\kappa \geq 1$ for the error-contaminated case, and also $\kappa = 0$ allowed in the error-free case. Accordingly, we impose the following assumption:
Comments on Assumptions (ref) and (ref) are in order. We assume that $\mathcal H_y = \mathbb{R}$ (resp.\ $\mathcal H_y = \mathcal H$) if $y_t$ is scalar-valued (resp.\ function-valued). In the function-valued case, $y_t$ may be defined not on $[a_1,a_2]$, as $x_t$ is, but on a different interval, say $[b_1,b_2]$; extending the subsequent theoretical results to this case is straightforward, and hence assuming $\mathcal H_y = \mathcal H$ entails no loss of generality. In Assumption (ref)(ref), we assume that \(d_N\) is finite. This condition has been widely employed in the literature on nonstationary functional time series and seems empirically relevant (Chang2016152,NSS,seoshang22). A wide class of functional time series satisfies this condition (see Remark (ref)). For convenience, we assume that \(d_N\) is known, even though it is not the case in most empirical applications. However, replacing \(d_N\) with a consistent estimator does not affect the asymptotic results to be developed. Moreover, we show in Section (ref) of the Supplementary Material that the variance-ratio testing procedure of NSS can be used in our setting, allowing for measurement errors. Assumption (ref)(ref) is employed to facilitate our theoretical analysis based on useful limit theorems in the existing literature (see, e.g., Bosq2000). Given Assumption (ref)(ref), Assumption (ref)(ref) does not appear restrictive, and some primitive sufficient conditions can be found in Bosq2000. Assumption (ref)(ref) is a technical condition required for our asymptotic analysis. Similar assumptions were employed by seo2020functional in the study of the functional principal component analysis (FPCA) for cointegrated functional time series. Sufficient, but not restrictive, conditions for the weak convergence results stated in Assumption (ref)(ref) can be found in, e.g., berkes2013weak and seo2020functional.
Assumption (ref) states requirements on the measurement errors. Although our primary focus is on the error-contaminated case (Assumption (ref)(ref)), we will also show how the theoretical results simplify in the error-free case (Assumption (ref)(ref)) when applying the standard covariance-based approach (see Remark (ref)). If \(e_t\) is serially independent and also independent of \(u_s^S\) and \(u_s^N\) for every \(s\) and \(t\), then, noting that (i) \(\mathbb{E}[e_t \otimes u_s^N] = \mathbb{E}[e_t \otimes u_s^S] = 0\) for all \(s\) and \(t\) in the considered scenario and (ii) \(e_t \otimes z_t\) is a Hilbert-valued random variable (see Theorems 2.7 and 2.16 of Bosq2000), the conditions in Assumption (ref)(ref) are satisfied under mild conditions. However, Assumption (ref)(ref) is not restricted to such a case and allows more general cases; specifically, \(e_t\) is assumed to be uncorrelated with \(z_{t+\ell}\) if \(|\ell|\) is sufficiently large, while no restriction is imposed on \(\mathbb{E}[e_t \otimes z_{t+\ell}]\) if \(|\ell|\) is small. That is, for our theoretical investigation of the proposed method, we only require each of \(e_t\), \(u_t^N\), and \(u_t^S\) to be uncorrelated with non-adjacent past or future measurement errors \(e_s\). This not only appears to be a mild assumption but is also reasonable for most empirical applications.
We introduce additional notation. When these quantities are well-defined, we let $\bm{\mathsf{\lambda}}_{j}[A]$ be the $j$-th largest eigenvalue of a compact operator $A$, $\bm{\mathsf{v}}_{j}[A] $ be the corresponding eigenvector, and $\bm{{\Pi}}_{j}[A] $ be the orthogonal projection onto $\operatorname{span}\{\bm{\mathsf{v}}_{j}[A]\}$; that is,
We also let $\Omega$ be defined by
Under Assumption (ref)(ref), the above is a well-defined bounded linear operator acting on $\mathcal H$ (see Section 2.3. of BSS2017). We hereafter let $\mathfrak F_t$ be the filtration given by
We first define the following operators for any nonnegative integer $\kappa\geq 0$:
Here, $\widehat{C}_{\kappa}$ is the so-called lag-$\kappa$ sample autocovariance operator and $\widehat{D}_{\kappa}$, by construction, is a nonnegative self-adjoint compact operator. As such, it allows the following spectral representation:
We then define its inverse on the restricted domain $\operatorname{ran} (\sum_{j=1}^{{\mathsf K}} \bm{{\Pi}}_{j}[\widehat{D}_{\kappa}])$ for ${\mathsf K} > 0$ as follows:
Our proposed estimator is constructed based on the following sample operator: for some random element $z_t$,
In the case where $\kappa=0$ (and thus $\widehat{C}_\kappa (\widehat{D}_\kappa)^{-1}_{\mathsf K}= \sum_{j=1}^{\mathsf K} \bm{\mathsf{\lambda}}^{-1}_{j}[\widehat{C}_{0}]\bm{{\Pi}}_{j}[\widehat{C}_{0}]$) and $z_t=y_t$, (ref) becomes identical to the standard FPCA-based estimator considered in the literature concerning stationary functional time series (see e.g., Park2012397). However, as the results of seong2021functional suggest, this estimator is affected by measurement errors and may not be a consistent estimator of $f$. In our context, $f$ on the subspace $\mathcal H^S =\operatorname{ran} {\rm P}^S$ is not generally consistently estimated, which results from the fact that the sample covariance ${\rm P}^S\widehat{C}_0{\rm P}^S$ suffers from non-negligible contamination by measurement errors (note that ${\rm P}^S\widehat{C}_0{\rm P}^S$ contains the component $T^{-1}\sum_{t=1}^T {\rm P}^S e_t\otimes {\rm P}^S e_t$, which is non-negligible) and thus is an inconsistent estimator of its true counterpart (i.e., $\mathbb{E}[{\rm P}^Sx_t\otimes {\rm P}^Sx_t]$). On the other hand, under Assumption (ref)(ref), we may deduce that ${\rm P}^S\widehat{C}_{\kappa}{\rm P}^S$ for $\kappa \geq 1$ does not suffer from such serious contamination. This is why we will mainly consider the case $\kappa \geq 1$ and construct our proposed estimator using the sample autocovariance operator. Computing $(\widehat{D}_\kappa)^{-1}_{\mathsf K}$ requires determining the number of retained eigenvectors, ${\mathsf K}$. Subsequently, we will require ${\mathsf K}$ to grow without bound depending on certain sample eigenvalues, but for now we assume only the following, required for the first few main results:
In our asymptotic analysis, we decompose $f$ as follows:
We then consistently estimate each summand. For our purposes, it is important to obtain consistent estimators of ${\rm P}^N$ and ${\rm P}^S$. We first show that such estimators can be obtained from the eigenvectors of $\widehat{D}_{\kappa}$. In the theorem below and hereafter, we let
In the case where there is no measurement error and $\kappa = 0$, we have $\bm{\mathsf{v}}_{j}[\widehat{D}_{0}] = \bm{\mathsf{v}}_{j}[\widehat{C}_{0}]$ and also $\Upsilon_T = 0$. This special case corresponds to Theorem 3.1 of seo2020functional, which concerns the FPCA of cointegrated functional time series, and Theorem (ref) can therefore be regarded as a suitable generalization of that result toward an autocovariance-based FPCA method that is robust to measurement errors. Theorem (ref) shows that the estimator $ \widehat{{\rm P}}_{\kappa}^{N}$ is super-consistent, and, as shown in our proof of Theorem (ref), the asymptotic bias remains and is of order $T^{-1}$; a similar result holds for $\widehat{{\rm P}}_{\kappa}^{S}$.
The projection estimators $ \widehat{{\rm P}}_{\kappa}^{N}$ and $ \widehat{{\rm P}}_{\kappa}^{S}$ give us a natural decomposition of $\widehat{D}_{\kappa}$. In the subsequent sections, we consider the decomposition of $\widehat{D}_{\kappa}$ in (ref) into the sum of $\widehat{D}_{\kappa}^{N}$ and $\widehat{D}_{\kappa}^{S}$ given (ref) below; this equation not only defines $\widehat{D}_{\kappa}^{N}$ and $\widehat{D}_{\kappa}^{S}$, but also highlights some of their key properties:
The properties above, along with the asymptotic properties of $\widehat{{\rm P}}_\kappa^{N}$ and $\widehat{{\rm P}}_\kappa^{S}$ in Theorem (ref), play a crucial role in the asymptotic analysis of our proposed estimator to be discussed.
Note that $f = f^N + f^S$, where $f^N$ captures how the persistent (nonstationary) component in $x_t$ affects $y_t$, while $f^S$ reflects the effect of the transitory (stationary) component. We propose an estimator for each of these two components, with the projections defined in (ref) playing a key role. Specifically, we propose estimators of $f$, $f^N$, and $f^S$, as follows:
where we note that
and these may be viewed as the inverses of $\widehat{D}_{\kappa}^{N}$ and $\widehat{D}_{\kappa}^{S}$ (see (ref)) in a restricted domain.
The following theorem establishes the consistency of $\widehat{f^N_{\kappa}}$ as an estimator of $f^N (= f{\rm P}^N)$ and details its limiting behavior: in the theorem below, $W^N$, $W^S$, and $W^u$ are defined as in Theorem (ref) and Assumption (ref), and recall that $\to_{p}$ denotes convergence in probability with respect to the usual operator norm for operator-valued sequences (see Section (ref)).
Theorem (ref) demonstrates that the proposed estimator $\widehat{f^N_{\kappa}}$ is a consistent estimator of $f^N$, and the asymptotic bias is of order $T^{-1}$; this result parallels that of the standard least squares-type estimator for the cointegrating relationship in the finite-dimensional case. We present remarks containing complementary results to Theorem (ref).
We next study the asymptotic properties of $\widehat{f^S_{\kappa}}$ as an estimator of $f^S(=f{\rm P}^S)$. As with standard FPCA-based estimators, our proposed estimator given in (ref) is defined on a finite-dimensional eigenspace of $\widehat{D}_\kappa^S$. Using the result that $\widehat{{\rm P}}^S_{\kappa} - {\rm P}^S = O_p(T^{-1})$ (see Theorem (ref)), we may deduce that $\widehat{D}_\kappa^S$ is a consistent estimator of $D_{\kappa}^S$, defined below:
Note that our estimator $\widehat{f^S_{\kappa}}$ is defined on a $({\mathsf K}-d_N)$-dimensional eigenspace of $\widehat{D}_{\kappa}$. For this estimator to be a consistent estimator of $f^S$ defined on the entire $\mathcal H^S$, we need some conditions on $D_{\kappa}^S$ and $f^S$. Moreover, it is also necessary to let ${\mathsf K}$ grow without bound. The required conditions are summarized below:
Assumption (ref)(ref) contains requirements similar to those employed by seong2021functional for functional linear models. The decision rule for ${\mathsf K}_S$ in Assumption (ref)(ref) adapts a commonly used approach, considered reasonable in practice for FPCA-based estimators (see, e.g., Section 3.1 and Remark 2 of the aforementioned paper). From the properties in (ref), we have
and hence computing ${\mathsf K}_{S}$ according to Assumption (ref)(ref) does not require additional calculation of eigenvalues associated with $\widehat{D}^S_{\kappa}$. We next give the asymptotic properties of $\widehat{f^S_{\kappa}}$ as an estimator of $f^S$. In the theorem below and hereafter, we let $\widetilde{C}_0^S = \mathbb{E}[{\rm P}^{S} \tilde x_{t} \otimes {\rm P}^{S} \tilde x_{t}]$, $\widetilde{C}_{u} = \mathbb{E}[\tilde{u}_t\otimes \tilde{u}_t]$, $\bm{{\tau}}_{j}[{D}_{\kappa}^S] = \max\{(\bm{\mathsf{\lambda}}_{j-1}[{D}_{\kappa}^S]-\bm{\mathsf{\lambda}}_{j}[{D}_{\kappa}^S])^{-1},(\bm{\mathsf{\lambda}}_{j}[{D}_{\kappa}^S]-\bm{\mathsf{\lambda}}_{j+1}[{D}_{\kappa}^S])^{-1}\}$,
where
Our next result studies the asymptotic properties of $\widehat{f^S_{\kappa}}$: in the theorem below, $N(0,A)$ denotes zero-mean Gaussian random element taking values in $\mathcal H_y$ with (co)variance $A$.
Even if the condition given by (ref) in Theorem (ref) requires that the eigenvalues of $D_{\kappa}^S$ are distinct, it does not place any other essential restrictions on the eigenstructure of $D_{\kappa}^S$. Given that $\sum_{j=1}^{{\mathsf K}_S}\bm{{\tau}}_{j}[D_{\kappa}^S]$ increases in ${\mathsf K}_S$ (and thus $\alpha^{-1})$, this condition merely requires $\alpha$ to decay to zero at a sufficiently slower rate. In fact, assumptions similar to (ref) are standard and widely used in the literature on functional linear models (see e.g., Park2012397,seong2021functional). Moreover, it is possible to relax the assumption of distinct eigenvalues in Theorem (ref) under a different set of assumptions, which is detailed in Remark (ref).
From Theorems (ref) and (ref), we know that the proposed estimator $\widehat{f}_{\kappa}$ is consistent under the employed assumptions, i.e., $\widehat{f}_{\kappa} \to_p f$. Moreover, as described by (ref), we find that our estimator $\widehat{f}_{\kappa}$ is asymptotically normal in a certain sense. However, unlike in a finite-dimensional setting, there are some limitations associated with the asymptotic normality given by (ref). First, $\widehat{f}_{\kappa}$ is centered at a random biased operator $f\widehat{{\rm P}}^{{\mathsf K}}_\kappa$, but not $f$, and (ii) the convergence is established in a pointwise manner at each point $\zeta \in \mathcal H$ but not uniformly over the entire space $\mathcal H$. As noted by seong2021functional, these limitations are in fact common in the literature concerning FPCA-based estimation of the functional linear model; see also Theorem 3.10 of Hu2016.
Obviously, (ref) may be used for inference on $f\widehat{{\rm P}}^{{\mathsf K}}_\kappa (\zeta)$, where $\widehat{{\rm P}}^{{\mathsf K}}_\kappa(\zeta)$ is naturally understood as the optimal approximation of $\zeta$ using the eigenvectors of $\widehat{D}_{\kappa}$. For example, when $\mathcal H_y = \mathcal H$ (and hence $y_t$ is function-valued), we may construct the 95% confidence interval of $\langle f\widehat{{\rm P}}^{{\mathsf K}}_\kappa (\zeta), \varphi \rangle$ for some $\varphi \in \mathcal H_y$ using the asymptotic normality result (ref), as follows::
where the unknown quantities ${\theta}_{{\mathsf K}_S}$ and $\widetilde{C}_{u}$ can be replaced by reasonable estimators that can be easily computed from our proposed estimator $\widehat{f}_{\kappa}$ without affecting asymptotic validity (see Corollary (ref) of the Supplementary Material). However, practitioners may want to avoid being interfered by a random projection $\widehat{{\rm P}}^{{\mathsf K}}_\kappa$ and implement a direct statistic inference on $\langle f (\zeta), \varphi\rangle$ rather than on $\langle f\widehat{{\rm P}}^{{\mathsf K}}_\kappa (\zeta), \varphi\rangle$. We will show that, under additional assumptions (requring “smoothness” of $f$ and $\zeta$), the asymptotic normality result (ref) still holds even when $f\widehat{{\rm P}}^{{\mathsf K}}_\kappa$ is replaced by $f$. This allows statistical inference without the influence of the random projection $\widehat{{\rm P}}^{{\mathsf K}}_\kappa$, as discussed in more detail in the next section.
We consider the following assumptions on ${\rm P}^Sx_t$, $f$ and $\zeta$, which are similar to the conditions employed by seong2021functional: below, for $j, \ell\geq 1$, $\varpi_{t}(j,\ell)=\langle {\rm P}_Sx_t, \bm{\mathsf{v}}_{j}[D_{\kappa}^S] \rangle\langle {\rm P}_Sx_{t-\kappa},\bm{\mathsf{v}}_{\ell}[E_{\kappa}^S] \rangle-\mathbb{E}[\langle {\rm P}_Sx_t,\bm{\mathsf{v}}_{j}[D_{\kappa}^S] \rangle\langle {\rm P}_Sx_{t-\kappa},\bm{\mathsf{v}}_{\ell}[E_{\kappa}^S]\rangle]$ and $E_{\kappa}^S = C^{S}_{\kappa}(C^S_{\kappa})^\ast$.
Assumption (ref)(ref) summarizes the technical conditions needed to establish the desired results. Similar requirements have been employed in the literature on functional linear models (see e.g., Hall2007,imaizumi2018,seong2021functional). Noting that, under this condition, $\bm{{\tau}}_{j}[D_{\kappa}^S]\leq c j^{\rho+1}$ and $\sum_{j=1}^{M} j^{\rho+1} = O(M^{\rho+2})$ for positive integer $M$, one may observe that, under Assumption (ref)(ref), the conditions given by (ref) may be replaced with the following sufficient condition: $T^{-1/2}\alpha^{-1/2} {\mathsf K}_S^{\rho+2} \to_p 0$. Given that $\alpha = a_1 T^{-a_2}$ for some $a_1 > 0$ and $a_2 \in (0, 1/2)$, Assumption (ref)(ref) requires that $\|f(\bm{\mathsf{v}}_{j}[D_{\kappa}^S])\|$ and $\langle \zeta, \bm{\mathsf{v}}_{j}[D_{\kappa}^S] \rangle$ decay to zero at a sufficiently fast rate as $j$ increases, implying that $f$ and $\zeta$ are sufficiently smooth with respect to the eigenvectors $\bm{\mathsf{v}}_{j}[D_{\kappa}^S]$.
Even if all the assumptions required for Theorem (ref) hold, $\widehat{f}_{\kappa}-f$ converges to a Gaussian random element at a rate depending on $\zeta$ and thus it is not generally possible to construct a uniform confidence band of $f$ from Theorem (ref) (see Remark (ref)). However, it may be possible to construct a local (or locally approximate) confidence band, which is naturally interpreted. We first note that $f(\zeta)$ may be understood as a partial effect on $y_t$ of a perturbation $\zeta$ in $x_t$, which is often of interest in practice. If $\mathcal H_y = \mathbb{R}$ and hence $y_t$ is real-valued, $f(\zeta)$ is a real-valued effect on $y_t$ of a perturbation $\zeta$, and in this case we may directly use (ref) for statistical inference by replacing $\widetilde{C}_{u}$ and $\theta_{{\mathsf K}_S}$ with their sample counterparts (see Corollary (ref) of the Supplementary Material). Now suppose that $\mathcal H_y = \mathcal H$. In this case, $f(\zeta)$ is a function defined on $[a_1,a_2]$, and we may construct a sequence of confidence intervals for local averages of $f(\zeta)$. Specifically, let $\mathcal I_j = (b_{j+1}- b_j)^{-1} 1\{u \in [b_j,b_{j+1}]\}$ for some $b_j$ and $b_{j+1}$ with $a_1 \leq b_j < b_{j+1} \leq a_2$. Then $\langle f(\zeta), \mathcal I_j \rangle = (b_{j+1}- b_j)^{-1} \int_{b_j}^{b_{j+1}} f(\zeta)(s) ds$ computes the local average of $f(\zeta)$ on the interval $[b_j,b_{j+1}]$. Using the results given in Theorem (ref), we know that
Note that $\langle\mathcal I_j, \widetilde{C}_{u} \mathcal I_j \rangle$ and $\theta_{{\mathsf K}_S}$ can also be replaced by their sample counterparts (Corollary (ref) of the Supplementary Material), enabling construction of an asymptotically valid confidence band for the local average via (ref). This can be applied to overlapping or non-overlapping sequences of intervals $\{\mathcal I_j\}_{j=1}^{M}$ with $\cup_{j=1}^M \mathcal I_j = [a_1,a_2]$. These confidence intervals are readily interpretable and can be useful in practice.
In the previous sections, we developed statistical inferential methods for the case where $\mathbb{E}[y_t] = \mathbb{E}[x_t] = 0$ for simplicity. However, in practice, a nonstationary time series may include a nonzero intercept, and hence our observations $x_t$ and $y_t$ may satisfy $\mathbb{E}[x_t] = \mu_x$ and $\mathbb{E}[y_t] = \mu_y$ for some unknown $\mu_x$ and $\mu_y$. To accommodate this scenario, one may consider the model with a deterministic term as follows:
where $\mu = f(\mu_x) - \mu_y \in \mathcal{H}$. With a straightforward modification, we can still achieve consistent estimation of $f$ and extend the statistical inference on $f(\zeta)$ given in Section (ref). More specifically, inference for this case can be implemented using the centered (demeaned) variables ${y}_{c,t}=y_t-\bar{y}_T$ and $\tilde{x}_{c,t}=x_t +e_t -\bar{x}_T-\bar{e}_T$, where $\bar{y}_T = T^{-1}\sum_{t=1}^T y_t$, and $\bar{x}_T$ and $\bar{e}_T$ are similarly defined. The proposed estimator is given as follows: $\widehat{f}_{c,\kappa} = \widehat{f}_{c,\kappa}^N + \widehat{f}_{c,\kappa}^S,$ where
where $\widehat{C}_{c,\kappa}$, $\widehat{D}_{c,\kappa}$, $\widehat{{\rm P}}_{c,\kappa}^N$, and $\widehat{{\rm P}}_{c,\kappa}^S$ are similarly computed as $\widehat{C}_{\kappa}$, $\widehat{D}_{\kappa}$, $\widehat{{\rm P}}_{\kappa}^N$, $\widehat{{\rm P}}_{\kappa}^S$, respectively, but with the centered variables. The consistency of the estimator can be established with only slight modifications, and the pointwise asymptotic normality can also be achieved as follows:
where $\theta_{c,{\mathsf K}_S}$ is defined similarly to $\theta_{{\mathsf K}_S}$ in (ref), but with $C_{\kappa}^S$ and the other operators ($\widetilde{C}_{0}^S$, $D_\kappa^S$, $(D_\kappa^S)_{{\mathsf K}_S}^{-1}$) depending on $C_{\kappa}^S$ being computed from the centered variables ${\rm P}^S\tilde{x}_{t} - \mathbb{E}[{\rm P}^S\tilde{x}_{t}]$ ($t=1,\ldots,T$) rather than ${\rm P}^S\tilde x_{t}$. As in the previous case, $\theta_{c,{\mathsf K}_{S}}$ and $\widetilde{C}_{u}$ can be replaced by their sample counterparts without affecting the asymptotic result in (ref), enabling statistical inference on $f(\zeta)$ in practice. A more detailed discussion including theoretical justification of these results is given in Section (ref) of the Supplementary Material.
We conducted simulation experiments to compare the autocovariance-based estimator with the standard covariance-based one. We compare estimator accuracy and coverage probabilities of confidence intervals for $\langle f(\zeta),\varphi \rangle$ for $\zeta \in \mathcal H$ and $\varphi \in \mathcal H_y$ across estimators. The simulation results support using the autocovariance-based approach with error-contaminated data. However, we postpone these details to Section (ref) of the Supplementary Material to focus on real-data analysis of the economic impact of climate change, which more effectively illustrates our methods' empirical relevance and usefulness.
We present an empirical study examining the economic impact of climate change using our proposed method. A large body of research in climate economics has shown that rising temperatures negatively affect economic growth (e.g., dell2012temperature; burke2015global; newell2021gdp; CruzHansberg2023), and estimating these effects is of central importance because they inform the economic costs of climate change and policy responses. Motivated by this literature, we apply our proposed method to global climate and economic data to demonstrate its empirical usefulness in providing statistical evidence for this relationship along with broader estimates of the economic impact of climate change.
We consider appropriate transformations of the probability densities of gross regional product (GRP) growth rates ($y_t$), as a measure of regional economic activity, and land temperature anomalies ($x_t$), commonly used as indicators of climate change. Because these time series are expected to be nonstationary, contaminated by measurement errors, and to exhibit nonzero unconditional means, we employ model (ref), with $u_t$ absorbing measurement errors in the dependent variable.
We use non-infilled gridded land temperature anomaly data from a collaborative product of the Climatic Research Unit at the University of East Anglia, the Met Office Hadley Centre, and the National Centre for Atmospheric Science (CRUTEM.5.0.2.0, osborn2021land). We estimate spatially distributed temperature anomaly densities for 1951–2019 using a Gaussian kernel with Silverman’s bandwidth. To avoid COVID-19-related distortions in $y_t$, data from 2020 onward are excluded. At each time $t$, the distribution’s support is restricted to the range containing 99% of the total probability mass, $[-5.80, 6.68]$, thereby excluding outliers as in chang2020evaluating.
For the GRP growth rates, we employ both the real GRP data of leonie_wenz_2023_7573249 for 1960–2019 and the real GDP in millions of 2021 international dollars, converted using Purchasing Power Parities, from the Conference Board Total Economy Database (TED) for 1950–2019.\footnote{The Conference Board Total Economy Database™ (April 2022) - Output, Labor and Labor Productivity, 1950-2022, downloaded from https://www.conference-board.org/data/economydatabase/total-economy-database-productivity on April 13, 2023.} leonie_wenz_2023_7573249 provide subnational economic output data for over 1,661 regions across 83 countries, enabling panel and cross-sectional regression analyses that reduce coverage bias and increase the number of observations. Based on this approach, we spatially disaggregate TED's country-level real GDP into the regional product levels from 1950 to 2019. Using these data, we estimate panel fixed-effects models to remove persistent regional heterogeneity and long-term structural changes, and then relate the residual component of regional growth to climate variables. Detailed procedures for spatial disaggregation, density generation (on support $[-0.105, 0.092]$, excluding extreme observations), and panel specifications are provided in Section (ref) of the Supplementary Material.
Figure (ref) shows the densities of land temperature anomalies and the temperature-related components of GRP growth rates, along with their first two central moments from 1951 to 2019. The mean of land temperature anomalies shows a persistent upward trend, with a rising standard deviation indicating greater variability. In contrast, the mean of temperature-related regional growth turned negative after the mid-1980s, while its standard deviation remained largely unchanged, suggesting stable dispersion in regional growth responses.
The literature notes that treating probability densities with support $[a_1,a_2]$ (without transformation) as standard Hilbert-valued elements is inadvisable (e.g., petersen2016), since densities do not form a linear space; this issue is particularly pronounced for nonstationary density-valued time series (seo2019cointegrated). To implement our framework, we apply the centered log-ratio (CLR) transformation of each density $g$, given by $g\mapsto \log g(u)-(a_2-a_1)^{-1} \int_{a_1}^{a_2} \log g(s)ds$ which, under regularity conditions, maps (in a bijective manner) the density $g$ into the subspace of $L^2[a_1,a_2]$ orthogonal to constant functions (egozcue2006hilbert), enabling direct application of our methods.\footnote{Since the CLR transformation involves $\log g(s)$, problems arise when $g(s)=0$. As is common in practice (e.g., seoshang22), this is avoided by adding a small constant to $g(s)$. In our study, the densities are constructed on a restricted domain excluding a few extreme values, so this issue does not occur.} The resulting CLR-transformed densities of GRP growth rates (resp.\ land temperature anomalies), generated from the raw data, are interpreted as measures of their true counterparts with measurement errors, and they are treated as functional data $y_t$ (resp. $\tilde{x}_t$) in (ref).\footnote{We assume that both time series include intercepts but no deterministic time trends, as in chang2020evaluating.}
We examine the nonstationarity of the CLR-transformed time series computed in the previous section and estimate the nonstationarity dimension $d_N$ of $\{x_t\}_{t\geq 1}$, an input to our inferential methods. We apply the variance-ratio testing procedure of NSS, which is shown to be robust to measurement errors in Section (ref) of the Supplementary Material. The procedure determines $d_N$ by testing $H_0: d_N = d_0$ against $H_1: d_N < d_0$ sequentially for $d_0 = d_{\max}, \ldots, 1$ until the null is not rejected for the first time, where $d_{\max}$ is a prespecified upper bound for $d_N$. The results with $d_{\max}=5$, presented in Table (ref), identify $d_N = 2$ at the standard 10% or 5% significance level. This testing procedure can also be useful for a different purpose in our empirical analysis. Our model requires that the stochastic trends of $y_t$ are explained by those of $x_t$ ($f^N$ captures this). A straightforward extension of this testing procedure can be used as a diagnostic check to see whether this is the case in practice, and the testing results are supportive of our empirical analysis (see Remark (ref) of the Supplementary Material for detailed discussion and results).
We present estimation results on the economic impact of climate change. Our main interest is in the slope parameter $f$ in the model with an intercept (ref). The estimation uses the CLR-transformed time series, with the threshold $\alpha$ set as in our simulation experiments (Section (ref) of the Supplementary Material), yielding ${\mathsf K}=4$ in this study. We use $d_N=2$, as estimated in Section (ref). Since $f^N$ can be viewed as capturing a stable long-run relationship between persistent stochastic trends, it is called the long-run response function, while $f^S$ is called the short-run response function, as in the literature (e.g., meierrieks2023temperature); accordingly, $f$ is referred to as the total-run response function.
We compute the estimators for $\kappa = 1$ and $\kappa = 0$ for comparison. The estimator $\widehat{f}_{\kappa}$ is an operator that maps a function to another. Although it can be visualized (since $f$ is Hilbert–Schmidt in the present setup, the estimated Hilbert–Schmidt kernel can be plotted in three dimensions), such a plot is unlikely to yield meaningful insights for practitioners concerned with the economic implications of climate-related scenarios or major events. Instead, we consider a functional change $\zeta$ in $x_t$, interpreted as a global warming shock to the world economy, and estimate its partial effect $f(\zeta)$ to quantify the economic damages resulting from the shock. The hypothetical global warming shock $\zeta$ can produce permanent effects, transitory effects, or both on regional economic growth. Permanent and transitory impacts are measured based on long-run (climate change) and short-run (interannual weather) variations of the functional changes, respectively. The estimated total-run response function thus illustrates how global warming collectively impacts the spatial distribution of regional growth rates. Note that while measurement errors do not affect the consistency of the long-run response function estimator $\widehat{f}_{\kappa}^N$ for either $\kappa=0$ or $\kappa=1$, they do affect the consistency of the short-run response function estimator $\widehat{f}_{\kappa}^S$ (see Remark (ref)). Thus, setting a positive $\kappa$ is necessary for robust statistical inference on the total-run response function.
To construct a representative global warming function, at first, we compute the mean difference between the first and second halves of the density estimates for the land temperature anomaly (hereafter, GW1). As shown in the left panel of Figure (ref), global warming can be conceptualized in statistical terms as a probabilistic shift from negative to positive anomalies, capturing the long-run distributional change in the Earth’s land surface temperature over the past 70 years. Previous studies estimate the break date for the northern hemisphere temperature anomaly at 1985 in the NASA dataset and 1984 in the HadCRUT3 dataset (estrada2013statistically; Estrada2019). Given the close similarity in statistical properties between the land temperature anomaly and the northern hemisphere series (chang2020evaluating), we adopt 1985 as a credible break date marking the onset of global warming in GW1. Of course, practitioners consider an alternative conceptualization of global warming. For example, one can define global warming as the distributional shift over the first and last five years of the sample period (hereafter, GW2); see the right panel of Figure (ref). In this setting, GW2 serves as a complementary measure, offering greater robustness to interannual variability and to uncertainties in the precise timing of the structural break. While Figure (ref) shows the probability densities, we use the model with CLR-transformed densities due to mathematical issues noted in the literature (e.g., petersen2016,seo2019cointegrated). Since the CLR map is bijective, we consider the CLR transformations of the densities in each panel of Figure (ref) and define $\zeta$ as the difference between the CLR-transformed densities. This is treated as a global warming shock in the model.
The estimate $\widehat{f}_{\kappa}(\zeta)$ captures the effect of the generated global warming shock on the CLR-transformed density of regional growth rates. Figure (ref) presents the estimated total-run $(\widehat{f}_{\kappa} (\zeta))$ and short-run $(\widehat{f}_{\kappa}^S (\zeta))$ responses to the considered global warming shock. Since the regressor is likely contaminated by measurement errors, statistical inference is conducted for the estimates with $\kappa = 1$, using the theoretical results in Theorem (ref) (and Corollary (ref) in the Supplementary Material). Specifically, the local confidence interval is obtained by estimating the pointwise standard error from the residual covariance within a one-grid bandwidth neighborhood and scaling it by the normal critical value at each point (see Section (ref) and (ref)). The 95% confidence intervals for the locally averaged response functions indicate that, while the short-run effects are statistically insignificant, global warming has a significant total-run impact on regional economic growth (potentially due to the limited sample size and the slower convergence rate of $\widehat{f}_{\kappa}^S$ compared to $\widehat{f}_{\kappa}^N$).
The downward slope of $\widehat{f}_{\kappa}(\zeta)$ indicates that global warming reduces the share of regions with high temperature-related economic growth while increasing the share with lower growth. In other words, as land temperatures rise, the distribution of regional growth shifts toward weaker outcomes. The slope is generally steeper under GW2 than under GW1, indicating that the magnitude of the climate-induced shift in regional growth outcomes is more pronounced when global warming is defined by the first-versus-last 5-year contrast. When measurement errors in the functional data are accounted for (with $\kappa = 1$), the slope of $\widehat{f}_{\kappa}(\zeta)$ becomes steeper at the right tail compared to the error-free case (with $\kappa = 0$). This discrepancy likely reflects bias from measurement errors, implying that the magnitude of climate-related economic impacts is underestimated when such errors are ignored.
From a practical perspective, it is more informative to visualize the distributional effect implied by $f(\zeta)$ in terms of changes in the probability density of GRP growth rates (noting that $f(\zeta)$ represents an effect on the CLR-transformed density). This is achieved by (i) fixing a reference density and its CLR transformation $y_{\text{ref}}$, and (ii) inverting the CLR-valued quantity $y_{\text{ref}} + \widehat{f}_{\kappa}(\zeta)$ back into the corresponding probability density, then comparing it with the reference density. For the inversion, the inverse CLR transformation $g(s) \mapsto \exp(g(s))/\int_{a_1}^{a_2} \exp(g(u))du$ is applied egozcue2006hilbert. Moreover, scaled global warming shocks $\zeta_{q} = q \zeta$ for $q \geq 0$ and their distributional effects are considered to examine how the reference density changes as the global warming shock intensifies or diminishes.
The left and middle panels of Figure (ref) show the result when the reference density is set to the average density of $y_t$ over the period 1951–1984, $q$ increases from $0$ to $1.5$, and $\zeta$ is constructed from GW1 or GW2. As $q$ increases, both shocks shift the mass of the distribution leftward and modestly widen it, reflecting lower average growth rates and greater dispersion across regions. The right panel summarizes these changes in terms of the first two moments. The mean declines approximately linearly with $q$, while the variance increases at an accelerating rate. Across all scales, GW2 produces more pronounced changes than GW1 in both the mean and variance, indicating a stronger impact on the central tendency and dispersion of regional growth rates. Taken together, these results suggest that stronger global-warming shocks are associated with slower average growth and increased dispersion, demonstrating the usefulness of our approach as a practical tool for policymakers to evaluate the adverse economic impacts of climate change.
This paper develops regression models for nonstationary and potentially error-contaminated functional time series and introduces a novel autocovariance-based inferential method. The methodology is broadly applicable to problems involving nonstationary functional data. Not only to illustrate our approach, but also for its intrinsic importance, we apply our methodology to assess the economic impact of climate change. Our analysis provides empirical evidence that global warming has a negative effect on regional economic growth.