Multiscale Comparison of Nonparametric Trend Curves

\heading{Multiscale Comparison}{of Nonparametric Trend Curves}

\authors{Marina Khismatullina\footnotemark[1]}{Erasmus University Rotterdam}{Michael Vogt\footnotemark[2]}{Ulm University} \footnotetext[1]{Corresponding author. Address: Erasmus School of Economics, Erasmus University Rotterdam, 3062 PA Rotterdam, Netherlands. Email: [email removed].} \footnotetext[2]{Address: Institute of Statistics, Department of Mathematics and Economics, Ulm University, 89081 Ulm, Germany. Email: [email removed].} \setcounter{footnote}{2}

Key words: Multiscale statistics; nonparametric regression; time series errors; shape constraints; strong approximations; anti-concentration bounds.

JEL classifications: C12; C14; C23; C38.

\numberwithin{equation}{section} \allowdisplaybreaks[1]

Introduction

The comparison of time trends is an important topic in time series analysis and has a wide range of applications in economics, finance and other fields such as climatology. In economics, one may wish to compare trends in real GDP Grier1989 or trends in long-term interest rates Christiansen1997 across different countries. In finance, it is of interest to compare the volatility trends of different stocks Nyblom2000. In climatology, researchers are interested in comparing the trending behaviour of temperature time series across different spatial locations KarolyWu2005.

In this paper, we develop new econometric methods for the comparison of time trends. Classically, time trends are modelled stochastically in econometrics; see e.g.\ Stock1988. Recently, however, there has been a growing interest in econometric models with deterministic time trends; see Cai2007, Atak2011, Robinson2012, ChenGaoLi2012, Zhang2012 and Hidalgo2014 among others. Following this recent development, we consider a general panel model with deterministic time trends. We observe a panel of $n$ time series $\mathcal{T}_i = \{ (Y_{it},\boldsymbol{X}_{it}): 1 \le t \le T \}$ for $1 \le i \le n$, where $Y_{it}$ are real-valued random variables and $\boldsymbol{X}_{it} = (X_{it,1},\ldots, X_{it, d})^\top$ are $d$-dimensional random vectors. Each time series $\mathcal{T}_i$ is modelled by the equation

for $1 \le t \le T$, where $\bm{\beta}_i$ is a $d \times 1$ vector of unknown parameters, $\boldsymbol{X}_{it}$ is a $d\times 1$ vector of individual covariates or controls, $m_i$ is an unknown nonparametric (deterministic) trend function defined on the rescaled time interval $[0,1]$, $\alpha_i$ is a fixed effect error term and $\mathcal{E}_i = \{ \varepsilon_{it}: 1 \le t \le T \}$ is a zero-mean stationary error process. A detailed description of model \tagform@{(ref)} together with the technical assumptions on the various model components can be found in Section (ref).

An important question in many applications is whether the nonparametric time trends $m_i$ in model \tagform@{(ref)} are all the same. This question can formally be addressed by a statistical test of the null hypothesis \[ H_0: m_1 = m_2 = \ldots = m_n \] against the alternative $H_1$: $m_i (u) \neq m_j(u)$ for some $u \in [0,1]$ and $i \ne j$. The problem of testing $H_0$ versus $H_1$ is well studied in the literature. In a restricted version of model \tagform@{(ref)} without covariates, it has been investigated in HaerdleMarron1990, Hall1990, DegrasWu2012 and ChenWu2018 among many others. Test procedures in a less restricted version of model \tagform@{(ref)} with covariates and homogeneous parameter vectors ($\beta_i = \beta$ for all $i$) have been derived, for instance, in Zhang2012 and Hidalgo2014.

Even though many different approaches to test $H_0$ versus $H_1$ have been developed over the years, most of them have a serious shortcoming: they are very uninformative. More specifically, they only allow to test whether the trend curves are the same or not. But they do not give any information on which curves are different and where (that is, in which time intervals) they differ. The following example illustrates the importance of this point: Suppose a group of researchers wants to test whether the GDP trend is the same in ten different countries. If the researchers use a test which is uninformative in the above sense and this test rejects $H_0$, they can only infer that there are at least two countries with a different trend. They do not obtain any information on which countries have a different GDP trend and where (that is, in which time periods) the trends differ. However, this is exactly the information that is relevant in practice. Without it, it is hardly possible to get a good economic understanding of the situation at hand. The main aim of the present paper is to construct a test procedure which provides the required information. More specifically, we construct a test which allows to make rigorous confidence statements about (i) which trend curves are different and (ii) \textit{where} (that is, in which time intervals) they differ.

Very roughly speaking, our approach works as follows: Let $\mathcal{I}_{u, h} := [u-h,u+h] \subseteq [0,1]$ be the rescaled time interval with midpoint $u \in [0,1]$ and length $2h$. We call $u$ the location and $h$ the scale of the interval $\mathcal{I}_{u, h}$. For any given location-scale pair $(u,h)$, let \[ H_0^{[i, j]}(u, h): m_i(w) = m_j(w) \text{ for all } w \in \mathcal{I}_{u, h} \] be the hypothesis that the two trend curves $m_i$ and $m_j$ are identical on the interval $\mathcal{I}_{u, h}$. Our procedure simultaneously tests $H_0^{[i, j]}(u, h)$ for a wide range of different location-scale pairs $(u,h)$ and all pairs $(i,j)$. As it takes into account multiple scales $h$, it is called a multiscale method. Our main theoretical result shows that under suitable technical conditions, the developed multiscale test controls the familywise error rate, that is, the probability of wrongly rejecting at least one null hypothesis $H_0^{[i, j]}(u, h)$. As we will see, this allows us to make simultaneous confidence statements of the following form for a given significance level $\alpha \in (0,1)$:

Hence, as desired, we can make confidence statements about which trend curves are different and where they differ.

To the best of our knowledge, there are no other test procedures available in the literature which allow to make simultaneous confidence statements of the form \tagform@{(ref)} in the context of the general model \tagform@{(ref)}. We are only aware of two other multiscale tests for the comparison of nonparametric time trends, both of which are restricted to a strongly simplified version of model \tagform@{(ref)} without covariates. The first one is a SiZer-type test by Park2009, for which theory has only been derived in the special case of $n=2$ time series. The second one is a multiscale test by KhismatullinaVogt2021, which was developed to detect differences between epidemic time trends in the context of the COVID-19 pandemic. Notably, our multiscale test can be regarded as an extension of the approach in KhismatullinaVogt2021. In the present paper, we go beyond the quite limited model setting in KhismatullinaVogt2021 and develop multiscale inference methods in the general model framework \tagform@{(ref)} which allows to deal with a wide range of economic and financial applications. Moreover, we develop a clustering algorithm based on our multiscale test which allows to detect groups of time series with the same trend.

Our multiscale test is constructed step by step in Section (ref), whereas its theoretical properties are laid out in Section (ref). The clustering algorithm is introduced and investigated in Section (ref). The proofs of all theoretical results are relegated to the technical Appendix. We complement the theoretical analysis of the paper by a simulation study and two application examples in Sections (ref) and (ref). In the first application example, we use our test and clustering methods to examine GDP growth data from different OECD contries. We in particular test the hypothesis that there is a common GDP trend in these countries and we cluster the countries into groups with the same trend. In the second example, we examine a long record of housing price data from different countries and compare the price trends in these countries by our methods.

The model framework

Notation

Throughout the paper, we adopt the following notation. For a vector $\mathbf{v} = (v_1, \ldots, v_m)\in\mathbb{R}^m$, we write $\|\mathbf{v}\|_q = \big(\sum_{i=1}^m v_i^q\big)^{1/q}$ to denote its $\ell_q$-norm and use the shorthand $\|\mathbf{v}\| = \|\mathbf{v}\|_2$ in the special case $q = 2$. For a random variable $V$, we define its $\mathcal{L}^q$-norm by $\|V\|_q = (\mathbb{E} |V|^q)^{1/q}$ and write $\|V\| := \|V\|_2$ in the case $q = 2$.

Let $\eta_t$ ($t \in \mathbb{Z}$) be independent and identically distributed ($\text{i.i.d.}$) random variables, write $\mathcal{F}_t = (\ldots, \eta_{t-1}, \eta_t)$ and let $g: \mathbb{R}^\infty \to \mathbb{R}$ be a measurable function such that $g(\mathcal{F}_t) = g(\ldots, \eta_{t-1}, \eta_t)$ is a properly defined random variable. Following Wu2005, we define the physical dependence measure of the process $\{g(\mathcal{F}_t)\}_{t=-\infty}^\infty$ by

where $\mathcal{F}_t^\prime = (\ldots, \eta_{-1}, \eta^\prime_0, \eta_1, \ldots, \eta_t)$ is a coupled version of $\mathcal{F}_t$ with $\eta_0^\prime$ being an i.i.d.\ copy of $\eta_0$. Evidently, $\delta_q(g, t)$ measures the dependency of the random variable $g(\mathcal{F}_t)$ on the innovation term $\eta_0$.

Model

We observe a panel of $n$ time series $\mathcal{T}_i = \{(Y_{it}, \boldsymbol{X}_{it}): 1 \le t \le T \}$ of length $T$ for $1 \le i \le n$. Each time series $\mathcal{T}_i$ satisfies the model equation

for $1 \le t \le T$, where $\bm{\beta}_i$ is a $d \times 1$ vector of unknown parameters, $\boldsymbol{X}_{it} = (X_{it,1},\ldots,X_{it,d})^\top$ is a $d\times 1$ vector of individual covariates, $m_i$ is an unknown nonparametric trend function defined on the unit interval $[0,1]$ with $\int_0^1 m_i(u) du = 0$ for all $i$, $\alpha_i$ is a (deterministic or random) intercept term and $\mathcal{E}_i = \{ \varepsilon_{it}: 1 \le t \le T \}$ is a zero-mean stationary error process. As common in nonparametric regression, the trend functions $m_i$ in model \tagform@{(ref)} depend on rescaled time $t/T$ rather than on real time $t$; see e.g.\ Robinson1989, Dahlhaus1997 and VogtLinton2014 for a discussion of rescaled time in nonparametric estimation. The condition $\int_0^1 m_i(u) du = 0$ is required for identification of the trend function $m_i$ in the presence of the intercept $\alpha_i$. We allow $\alpha_i$ to be correlated with the covariates $\boldsymbol{X}_{it}$ in an arbitrary way. Hence, $\alpha_i$ can be regarded as a fixed effect error term. We do not impose any restrictions on the dependence between the fixed effects $\alpha_i$ across $i$. Similarly, the covariates $\boldsymbol{X}_{it}$ are allowed to be dependent across $i$ in an arbitrary way. As a consequence, the $n$ time series $\mathcal{T}_i$ in our panel can be correlated with each other in various ways. In contrast to $\alpha_i$ and $\boldsymbol{X}_{it}$, the error terms $\varepsilon_{it}$ are assumed to be independent across $i$. Technical conditions regarding the model are discussed below. Throughout the paper, we restrict attention to the case where the number of time series $n$ in model \tagform@{(ref)} is fixed. Extending our theoretical results to the case where $n$ grows with the time series length $T$ is a possible topic for further research.

Assumptions

The error processes $\mathcal{E}_i = \{ \varepsilon_{it}: 1 \le t \le T\}$ satisfy the following conditions.

Assumption (ref) implies that the error processes $\mathcal{E}_i$ are stationary and causal (in the sense that $\varepsilon_{it}$ does not depend on future innovations $\eta_{is}$ with $s>t$). The class of error processes that satisfy condition (ref) is very large. It includes linear processes, nonlinear transformations thereof, as well as a large variety of nonlinear processes such as Markov chain models and nonlinear autoregressive models Wu2016. Following Wu2005, we impose conditions on the dependence structure of the error processes $\mathcal{E}_i$ in terms of the physical dependence measure $\delta_q(g_i, t)$ defined in \tagform@{(ref)}. In particular, we assume the following:

For fixed $i$ and $t$, the expression $\sum\nolimits_{s \ge t} \delta_q(g_i, s)$ measures the cumulative effect of the innovation $\eta_0$ on the variables $\varepsilon_{it}, \varepsilon_{it+1},\ldots$ in terms of the $\mathcal{L}^q$-norm. Condition (ref) puts some restrictions on the decay of $\sum\nolimits_{s \ge t} \delta_q(g_i, s)$ (as a function of $t$). It is fulfilled by a wide range of stationary processes $\mathcal{E}_i$. For a detailed discussion of (ref)--(ref) and some examples of error processes that satisfy these conditions, see KhismatullinaVogt2020.

The covariates $\boldsymbol{X}_{it} = (X_{it,1},\ldots, X_{it, d})^\top$ are assumed to have the following properties.

Assumption (ref) guarantees that the process $\{ \boldsymbol{X}_{it}: 1 \le t \le T \}$ is stationary and causal for each $i$. Assumption (ref) restricts the serial dependence of the process $\{ \boldsymbol{X}_{it}: 1 \le t \le T \}$ for each $i$ in terms of the physical dependence measure.

We finally impose some assumptions on the relationship between the covariates and the errors and on the trend functions $m_i$.

The multiscale test

In this section, we develop a multiscale test for the comparison of the trend curves $m_i$ in model \tagform@{(ref)}. More specifically, we construct a multiscale test of the null hypothesis \[ H_0: m_1 = m_2 = \ldots = m_n. \] The test is designed to be as informative as possible: It does not only allow to say whether $H_0$ is violated. It also gives information on which types of violation occur. In particular, it allows to infer which trends are different and where (that is, in which time intervals) they differ.

Our strategy to construct the test is as follows: For any given location-scale point $(u,h)$, let \[ H_0^{[i, j]}(u, h): m_i(w) = m_j(w) \text{ for all } w \in \mathcal{I}_{u, h} \] be the hypothesis that $m_i$ and $m_j$ are identical on the interval $\mathcal{I}_{u, h} := [u-h,u+h] \subseteq [0,1]$. $H_0^{[i, j]}(u, h)$ can be viewed as a local null hypothesis that characterises the behaviour of two trend functions locally on the interval $\mathcal{I}_{u, h}$, whereas $H_0$ is the global null hypothesis that is concerned with the comparison of all trends on the whole unit interval $[0, 1]$. We consider a large family of time intervals $\mathcal{I}_{u, h}$ that fully cover the unit interval. More formally, we consider all intervals $\mathcal{I}_{u, h}$ with $(u,h) \in \mathcal{G}_T$, where $\mathcal{G}_T$ is a set of location-scale points $(u,h)$ with the property that $\bigcup_{(u, h) \in \mathcal{G}_T} \mathcal{I}_{u, h} = [0,1]$. Technical conditions on the set $\mathcal{G}_T$ are given below. The global null $H_0$ can now be formulated as

This formulation shows that we can test $H_0$ by a procedure which simultaneously tests the local hypothesis $H_0^{[i, j]}(u, h)$ for all $1 \le i < j \le n$ and all $(u,h) \in \mathcal{G}_T$.

In what follows, we design such a simultaneous test procedure step by step. To start with, we introduce some auxiliary estimators in Section (ref). We then construct the test statistics in Section (ref) and set up the test procedure in Section (ref). Finally, Section (ref) explains how to implement the procedure in practice.

Preliminary steps

If the fixed effects $\alpha_i$ and the coefficient vectors $\bm{\beta}_i$ were known, our testing problem would be greatly simplified. In particular, we could consider the model

which is a standard nonparametric regression equation. As $\alpha_i$ and $\bm{\beta}_i$ are not observed in practice, we construct estimators $\widehat{\alpha}_i$ and $\widehat{\bm{\beta}}_i$ and replace the unknown variables $Y_{it}^\circ$ by the approximations

In what follows, we explain how to construct estimators of $\alpha_i$ and $\bm{\beta}_i$.

In order to estimate the parameter vector $\bm{\beta}_i$ for a given $i$, we consider the time series $\Delta \mathcal{T}_i = \{(\Delta Y_{it}, \Delta \boldsymbol{X}_{it}): 2 \leq t \leq T\}$ of the first differences $\Delta Y_{it} = Y_{it} - Y_{i t-1}$ and $\Delta \boldsymbol{X}_{it} = \boldsymbol{X}_{it} - \boldsymbol{X}_{it-1}$. We can write

where $ \Delta \varepsilon_{it} = \varepsilon_{it} - \varepsilon_{i t-1}$. Since $m_i$ is Lipschitz by Assumption (ref), we can use the fact that $ |m_i ( \frac{t}{T} ) - m_i (\frac{t-1}{T}) | = O(\frac{1}{T})$ to get that

This suggests to estimate $\bm{\beta}_i$ by applying least squares methods to equation \tagform@{(ref)}, treating $\Delta \boldsymbol{X}_{it}$ as the regressors and $\Delta Y_{it}$ as the response variable. As a result, we obtain the least squares estimator

In Lemma (ref) in the Appendix, we show that $\widehat{\bm{\beta}}_i - \bm{\beta}_i = O_p(T^{-1/2})$ under our assumptions. Given $\widehat{\bm{\beta}}_i$, we next estimate the fixed effect term $\alpha_i$ by

This is a reasonable estimator of $\alpha_i$ for the following reason: For each $i$, it holds that

since $\frac{1}{T}\sum_{i=1}^T \varepsilon_{it} = O_p(T^{-1/2})$ by the law of large numbers, $\frac{1}{T}\sum_{i=1}^T m_i(t/T) = O(T^{-1})$ due to Lipschitz continuity of $m_i$ and the normalisation $\int_{0}^1 m_i(u)du = 0$, $\frac{1}{T}\sum_{t=1}^T \boldsymbol{X}_{it} = O_p(1)$ by Chebyshev's inequality and $\widehat{\bm{\beta}}_i - \bm{\beta}_i = O_p (T^{-1/2})$.

In order to construct our multiscale test, we do not only require the variables $\widehat{Y}_{it}$ and thus the estimators $\widehat{\bm{\beta}}_i$ and $\widehat{\alpha}_i$. We also need an estimator of the long-run error variance $\sigma_i^2 = \sum\nolimits_{\ell=-\infty}^{\infty} \textnormal{Cov}(\varepsilon_{i0}, \varepsilon_{i\ell})$ for each $i$. Throughout the paper, we assume that the long-run variance $\sigma_i^2$ does not depend on $i$, that is $\sigma_i^2 = \sigma^2$ for all $i$. For technical reasons, we nevertheless use a different estimator of $\sigma_i^2 = \sigma^2$ for each $i$. Let $\widehat{\sigma}_i^2$ be an estimator of $\sigma_i^2$ which is computed from the $i$-th time series $\{ \widehat{Y}_{it}: 1 \le t \le T \}$, that is, let $\widehat{\sigma}_i^2 = \widehat{\sigma}_i^2(\widehat{Y}_{i1},\ldots,\widehat{Y}_{iT})$ be a function of the variables $\widehat{Y}_{it}$ for $1 \le t \le T$. Our theory works with any estimator $\widehat{\sigma}_i^2$ which has the property that $\widehat{\sigma}_i^2 = \sigma_i^2 + o_p(\rho_T)$ for each $i$, where $\rho_T$ slowly converges to zero, in particular, $\rho_T = o(1/\log T)$. We now briefly discuss two possible choices of $\widehat{\sigma}_i^2$.

Following Kim2016, we can estimate $\sigma_i^2$ by a variant of the subseries variance estimator, which was first proposed by Carlstein1986 and then extended by WuZhao2007. Formally, we set

where $s_T$ is the length of the subseries and $M = \lfloor T/s_T\rfloor$ is the largest integer not exceeding $T/s_T$. As per the optimality result in Carlstein1986, we set $s_T \asymp T^{1/3}$. When implementing $\widehat{\sigma}_i^2$, we in particular choose $s_T = \lfloor T^{1/3}\rfloor$. According to Lemma (ref) in the Appendix, $\widehat{\sigma}_i^2$ is an asymptotically consistent estimator of $\sigma_i^2$ with the rate of convergence $O_p(T^{-1/3})$.

The subseries estimator $\widehat{\sigma}_i^2$ is completely nonparametric: it does not presuppose a specific time series model for the error process $\mathcal{E}_i = \{ \varepsilon_{it}: 1 \le t \le T\}$ but merely requires $\mathcal{E}_i$ to fulfill general weak dependence conditions. In practice, however, it often makes sense to impose additional structure on the error process $\mathcal{E}_i$. A very popular yet general error model is the class of $\text{AR}(\infty)$ processes. A difference-based estimator of the long-run error variance for this error model has been developed in KhismatullinaVogt2020 and can be easily adapted to the situation at hand. A detailed discussion of the estimator and a comparison with other long-run variance estimators can be found in Section 4 of KhismatullinaVogt2020.

Construction of the test statistics

A statistic to test the local null hypothesis $H_0^{[i, j]}(u, h)$ for a given location-scale point $(u,h)$ and a given pair of time series $(i,j)$ can be constructed as follows: Let

be a weighted average of the differenced variables $\widehat{Y}_{it} - \widehat{Y}_{jt}$. In particular, let $w_{t,T}(u, h)$ be local linear kernel weights of the form

where \[ \Lambda_{t,T}(u, h) = K\Big(\frac{\frac{t}{T}-u}{h}\Big) \Big[ S_{T,2}(u, h) - \Big(\frac{\frac{t}{T}-u}{h}\Big) S_{T,1}(u, h) \Big], \] $S_{T,\ell}(u, h) = (Th)^{-1} \sum\nolimits_{t=1}^T K(\frac{\frac{t}{T}-u}{h}) (\frac{\frac{t}{T}-u}{h})^\ell$ for $\ell = 1,2$ and $K$ is a kernel function with the following property:

The kernel average $\widehat{\psi}_{ij, T}(u, h)$ can be regarded as a measure of distance between the two trend curves $m_i$ and $m_j$ on the interval $\mathcal{I}_{u, h} = [u-h, u+h]$. We may thus use a normalized (absolute) version of the kernel average $\widehat{\psi}_{ij, T}(u, h)$, in particular, the term

as a statistic to test the local null $H_0^{[i, j]}(u, h)$. Our aim is to test $H_0^{[i, j]}(u, h)$ simultaneously for a wide range of location-scale points $(u,h)$ and all possible pairs of time series $(i,j)$. To take into account that we are faced with a simultaneous test problem, we replace the statistics from \tagform@{(ref)} by additively corrected versions of the form

where $\lambda(h) = \sqrt{2 \log \{ 1/(2h) \}}$. The scale-dependent correction term $\lambda(h)$ was first introduced in the multiscale approach of DuembgenSpokoiny2001 and has been used in various contexts since then. We refer to KhismatullinaVogt2020 for a detailed discussion of the idea behind this additive correction.

In order to test the global null hypothesis $H_0$, we aggregate the test statistics $\widehat{\psi}^0_{ij, T}(u, h)$ by taking their maximum over all location-scale points $(u,h) \in \mathcal{G}_T$ and all pairs of time series $(i,j)$ with $1 \le i < j \le n$. This leads to the multiscale test statistic \[ \widehat{\Psi}_{n,T} = \max_{(u,h) \in \mathcal{G}_T, 1 \le i < j \le n} \widehat{\psi}^0_{ij, T}(u, h). \] For our theoretical results, we suppose that the set of location-scale points $\mathcal{G}_T$ is a subset of $\mathcal{G}_T^{\text{full}} = \{ (u, h): \mathcal{I}_{u, h} = [u-h,u+h] \subseteq [0,1]$ with $u = t/T$ and $h = s/T$ for some $1 \le t, s \le T$ and $h \in [h_{\min},h_{\max}] \}$ which fulfills the following conditions:

Assumptions (ref) and (ref) place relatively mild restrictions on the set $\mathcal{G}_T$: (ref) allows $\mathcal{G}_T$ to be very large compared to the sample size $T$. In particular, $\mathcal{G}_T$ may grow as any polynomial of $T$. (ref) allows $\mathcal{G}_T$ to contain intervals $[u-h,u+h]$ of many different scales $h$, ranging from very small intervals of scale $h_{\min}$ to substantially larger intervals of scale $h_{\max}$.

The test procedure

Let $\mathcal{M} = \{ (u,h,i,j): (u,h) \in \mathcal{G}_T \text{ and } 1 \le i < j \le n \}$ be the collection of all location-scale points $(u,h)$ and all pairs of time series $(i,j)$ under consideration. Moreover, let $\mathcal{M}_0 \subseteq \mathcal{M}$ be the set of tuples $(u,h,i,j)$ for which $H_0^{[i, j]}(u, h)$ is true. For a given significance level $\alpha \in (0,1)$, our multiscale test is carried out as follows:

We now construct a critical value $q_{n,T}(\alpha)$ which controls the FWER at level $\alpha$. Let $q_{n,T}(\alpha)$ be the $(1-\alpha)$-quantile of the multiscale statistic $\widehat{\Psi}_{n,T}$ under the global null $H_0: m_1 = m_2 = \ldots = m_n$. Since

with $\mathbb{P}_{H_0}$ denoting the probability under $H_0$, this choice of $q_{n,T}(\alpha)$ indeed controls the FWER at the desired level. However, this choice is not computable in practice. We thus replace it by a suitable approximation: Suppose we could perfectly estimate $\bm{\beta}_i$ and $\sigma_i^2$, that is, $\widehat{\bm{\beta}}_i = \bm{\beta}_i$ and $\widehat{\sigma}_i^2 = \sigma_i^2$ for all $i$. Then under $H_0$, the test statistics $\widehat{\psi}^0_{ij, T}(u, h)$ would simplify to

where $\bar{\varepsilon}_i = \bar{\varepsilon}_{i,T} := T^{-1} \sum_{t=1}^T \varepsilon_{it}$ denotes the empirical average of the variables $\varepsilon_{i1},\ldots,\varepsilon_{iT}$. Replacing the error terms $\varepsilon_{it}$ in \tagform@{(ref)} by normally distributed variables leads to the statistics

with

where $Z_{it}$ are independent standard normal random variables for $1 \le t \le T$ and $1 \le i \le n$ and, as before, we use the shorthand $\bar{Z}_i = \bar{Z}_{i,T} := T^{-1} \sum_{t=1}^T Z_{it}$. With this notation, we define

which can be regarded as a Gaussian version of the test statistic $\widehat{\Psi}_{n,T}$ under the null $H_0$ (in the idealized case with $\widehat{\bm{\beta}}_i = \bm{\beta}_i$ and $\widehat{\sigma}_i^2 = \sigma_i^2$ for all $i$). We now choose $q_{n,T}(\alpha)$ to be the $(1-\alpha)$-quantile of $\Phi_{n,T}$.

Implementation of the test in practice

In practice, we implement the test procedure as follows for a given significance level $\alpha \in (0, 1)$:

Theoretical properties of the multiscale test

To start with, we investigate the auxiliary statistic

where

and

with $\bar{\varepsilon}_i = \bar{\varepsilon}_{i,T} := T^{-1} \sum_{t=1}^T \varepsilon_{it}$ and $\bar{\boldsymbol{X}}_{i} = \bar{\boldsymbol{X}}_{i, T} := T^{-1}\sum_{t=1}^T \boldsymbol{X}_{it}$. By construction, it holds that $\widehat{\phi}_{ij, T}(u, h) = \widehat{\psi}_{ij, T}(u, h)$ under $H_0^{[i, j]}(u, h)$. Hence, the auxiliary statistic $\widehat{\Phi}_{n,T}$ is identical to the multiscale test statistic $\widehat{\Psi}_{n,T}$ under the global null $H_0$. Consequently, in order to determine the distribution of the test statistic $\widehat{\Psi}_{n,T}$ under $H_0$, it suffices to study the auxiliary statistic $\widehat{\Phi}_{n,T}$. The following theorem shows that the distribution of $\widehat{\Phi}_{n,T}$ is close to the distribution of the Gaussian statistic $\Phi_{n,T}$ introduced in \tagform@{(ref)}.

Theorem (ref) is key for deriving theoretical properties of our multiscale test. Its proof is provided in the Appendix.

With the help of Theorem (ref), we now examine the theoretical properties of our multiscale test developed in Section (ref). The following proposition shows that our test of the global null $H_0$ has correct (asymptotic) size.

The next proposition characterises the power of the test against a certain class of local alternatives. To formulate it, we consider a sequence of pairs of functions $m_ i := m_{i,T}$ and $m_ j := m_{j,T}$ that depend on the time series length $T$ and that are locally sufficiently far from each other.

We now turn attention to the local null hypotheses $H_0^{[i, j]}(u, h)$. As already defined above, let $\mathcal{M} = \{ (u,h,i,j): (u,h) \in \mathcal{G}_T \text{ and } 1 \le i < j \le n \}$ be the collection of location-scale points $(u,h)$ and pairs of time series $(i,j)$ under consideration. Moreover, let $\mathcal{M}_0 \subseteq \mathcal{M}$ be the set of tuples $(u,h,i,j)$ for which $H_0^{[i, j]}(u, h)$ is true. Since we test $H_0^{[i, j]}(u, h)$ simultaneously for all $(u,h,i,j) \in \mathcal{M}$, we would like to bound the probability of making at least one false discovery. More formally, we would like to control the family-wise error rate

at level $\alpha$, where $\mathcal{S}^{[i, j]}(\alpha)$ is the set of intervals $\mathcal{I}_{u, h}$ for which the null hypothesis $H_0^{[i, j]}(u, h)$ is rejected. The following result shows that our test procedure asymptotically controls the FWER at level $\alpha$.

The statement of Proposition (ref) can obviously be reformulated as follows:

Consequently, we can make the following simultaneous confidence statement:

Put differently:

Clustering

Consider a situation in which the null hypothesis $H_0: m_1 = m_2 = \ldots = m_n$ is violated. Even though some of the trend functions are different in this case, part of them may still be the same. Put differently, there may be groups of time series which have the same time trend. Formally speaking, we define a group structure as follows: There exist sets or groups of time series $G_1,\ldots, G_N$ with $N \le n$ and $\{1,\ldots, n\} = \mathbin{\dot{\bigcup}}_{\ell=1}^{N} G_\ell$ such that for each $1 \le \ell \le N$, \[ m_i = f_\ell \quad \text{for all } i \in G_\ell, \] where $f_\ell$ are group-specific trend functions. Hence, the time series which belong to the group $G_\ell$ all have the same time trend $f_\ell$. Throughout the section, we suppose that the group-specific trend functions $f_\ell$ have the following properties:

In many applications, it is natural to suppose that there is a group structure in the data. In this case, a particular interest lies in estimating the unknown groups from the data at hand. In what follows, we combine our multiscale methods with a clustering algorithm to achieve this. More specifically, we use the multiscale statistics $\max_{(u, h) \in \mathcal{G}_T}\hat{\psi}^0_{ij, T}(u, h)$ calculated for each $i$ and $j$ as distance measures which are fed into a hierarchical clustering algorithm.

To describe the algorithm, we first need to introduce the notion of a dissimilarity measure: Let $S \subseteq \{1,\ldots, n\}$ and $S^\prime \subseteq \{1,\ldots, n\}$ be two sets of time series from our sample. We define a dissimilarity measure between $S$ and $S^\prime$ by setting

This is commonly called a complete linkage measure of dissimilarity. Alternatively, we may work with an average or a single linkage measure. We now combine the dissimilarity measure $\widehat{\Delta}$ with a hierarchical agglomerative clustering (HAC) algorithm which proceeds as follows:

Step $0$ (Initialisation): Let $\widehat{G}_i^{[0]} = \{ i \}$ denote the $i$-th singleton cluster for $1 \le i \le n$ and define $\{\widehat{G}_1^{[0]},\ldots,\widehat{G}_n^{[0]} \}$ to be the initial partition of time series into clusters.

Step $r$ (Iteration): Let $\widehat{G}_1^{[r-1]},\ldots,\widehat{G}_{n-(r-1)}^{[r-1]}$ be the $n-(r-1)$ clusters from the previous step. Determine the pair of clusters $\widehat{G}_{\ell}^{[r-1]}$ and $\widehat{G}_{{\ell}^\prime}^{[r-1]}$ for which

\[ \widehat{\Delta}(\widehat{G}_{\ell}^{[r-1]},\widehat{G}_{{\ell}^\prime}^{[r-1]}) = \min_{1 \le k < k^\prime \le n-(r-1)} \widehat{\Delta}(\widehat{G}_{k}^{[r-1]},\widehat{G}_{k^\prime}^{[r-1]}) \] and merge them into a new cluster.

Iterating this procedure for $r = 1,\ldots, n-1$ yields a tree of nested partitions $\{\widehat{G}_1^{[r]},\ldots$ $\ldots,\widehat{G}_{n-r}^{[r]}\}$, which can be graphically represented by a dendrogram. Roughly speaking, the HAC algorithm merges the $n$ singleton clusters $\widehat{G}_i^{[0]} = \{ i \}$ step by step until we end up with the cluster $\{1,\ldots, n\}$. In each step of the algorithm, the closest two clusters are merged, where the distance between clusters is measured in terms of the dissimilarity $\widehat{\Delta}$. We refer the reader to Section 14.3.12 in HastieTibshiraniFriedman2009 for an overview of hierarchical clustering methods.

When the number of groups $N$ is known, we estimate the group structure $\{G_1,\ldots, G_N\}$ by the $N$-partition $\{\widehat{G}_1^{[n-N]},\ldots,\widehat{G}_{N}^{[n-N]}\}$ produced by the HAC algorithm. When $N$ is unknown, we estimate it by the $\widehat{N}$-partition $\{\widehat{G}_1^{[n-\widehat{N}]},\ldots,\widehat{G}_{\widehat{N}}^{[n-\widehat{N}]}\}$, where $\widehat{N}$ is an estimator of $N$. The latter is defined as \[ \widehat{N} = \min \Big\{ r = 1,2,\ldots \Big| \max_{1 \le \ell \le r} \widehat{\Delta} \big( \widehat{G}_\ell^{[n-r]} \big) \le q_{n,T}(\alpha) \Big\}, \] where we write $\widehat{\Delta}(S) = \widehat{\Delta}(S,S)$ for short and $q_{n,T}(\alpha)$ is the $(1-\alpha)$-quantile of $\Phi_{n,T}$ defined in Section (ref).

The following proposition summarises the theoretical properties of the estimators $\widehat{N}$ and $\{ \widehat{G}_1,\ldots,\widehat{G}_{\widehat{N}} \}$, where we use the shorthand $\widehat{G}_\ell = \widehat{G}_\ell^{[n-\widehat{N}]}$ for $1 \le \ell \le \widehat{N}$.

This result allows us to make statistical confidence statements about the estimated clusters $\{ \widehat{G}_1,\ldots,\widehat{G}_{\widehat{N}} \}$ and their number $\widehat{N}$. In particular, we can claim with asymptotic confidence at least $1 - \alpha$ that the estimated group structure is identical to the true group structure. The proof of Proposition (ref) can be found in the Appendix.

Our multiscale methods do not only allow us to compute estimators of the unknown groups $G_1,\ldots, G_N$ and their number $N$. They also provide information on the locations where two group-specific trend functions $f_\ell$ and $f_{\ell^\prime}$ differ from each other. To turn this claim into a mathematically precise statement, we need to introduce some notation. First of all, note that the indexing of the estimators $\widehat{G}_1,\ldots,\widehat{G}_{\widehat{N}}$ is completely arbitrary. We could, for example, change the indexing according to the rule $\ell \mapsto \widehat{N} - \ell + 1$. In what follows, we suppose that the estimated groups are indexed such that $P( \widehat{G}_\ell = G_\ell \text{ for all } \ell ) \ge (1-\alpha) + o(1)$. Proposition (ref) implies that this is possible without loss of generality. Keeping this convention in mind, we define the sets \[ \mathcal{A}_{n,T}^{[\ell,\ell^\prime]}(\alpha)= \Big\{ (u, h) \in \mathcal{G}_T: \widehat{\psi}^0_{ij, T}(u, h) > q_{n,T}(\alpha) \text{ for some } i \in \widehat{G}_\ell, j \in \widehat{G}_{\ell^\prime} \Big\} \] and \[ \mathcal{S}^{[\ell,\ell^\prime]}_{n, T}(\alpha) = \big\{ \mathcal{I}_{u, h} = [u-h, u+h]: (u, h) \in \mathcal{A}_{n,T}^{[\ell,\ell^\prime]}(\alpha) \big\} \] for $1 \le \ell < \ell^\prime \le \widehat{N}$. An interval $\mathcal{I}_{u, h}$ is contained in $\mathcal{S}^{[\ell,\ell^\prime]}_{n, T}(\alpha)$ if our multiscale test indicates a significant difference between the trends $m_i$ and $m_j$ on the interval $\mathcal{I}_{u, h}$ for some $i \in \widehat{G}_\ell$ and $j \in \widehat{G}_{\ell^\prime}$. Put differently, $\mathcal{I}_{u, h} \in \mathcal{S}^{[\ell,\ell^\prime]}_{n, T}(\alpha)$ if the test suggests a significant difference between the trends of the $\ell$-th and the $\ell^\prime$-th group on the interval $\mathcal{I}_{u, h}$. We further let \[ E_{n,T}^{[\ell,\ell^\prime]}(\alpha) = \Big\{ \forall \mathcal{I}_{u, h} \in \mathcal{S}^{[\ell,\ell^\prime]}_{n, T}(\alpha): f_\ell(v) \ne f_{\ell^\prime}(v) \text{ for some } v \in \mathcal{I}_{u, h} = [u-h, u+h] \Big\} \] be the event that the group-specific time trends $f_\ell$ and $f_{\ell^\prime}$ differ on all intervals $\mathcal{I}_{u, h} \in \mathcal{S}^{[\ell,\ell^\prime]}_{n, T}(\alpha)$. With this notation at hand, we can make the following formal statement whose proof is given in the Appendix.

According to Proposition (ref), the sets $\mathcal{S}^{[\ell,\ell^\prime]}_{n, T}(\alpha)$ allow us to locate, with a pre-specified confidence, time intervals where the group-specific trend functions $f_\ell$ and $f_{\ell^\prime}$ differ from each other. In particular, we can claim with asymptotic confidence at least $1 - \alpha$ that the trend functions $f_\ell$ and $f_{\ell^\prime}$ differ on all intervals $\mathcal{I}_{u, h} \in \mathcal{S}^{[\ell,\ell^\prime]}_{n, T}(\alpha)$.\footnote{This statement remains to hold true when the sets of intervals $\mathcal{S}^{[\ell,\ell^\prime]}_{n, T}(\alpha)$ are replaced by the corresponding sets of minimal intervals.}

Simulations

In this section, we investigate our testing and clustering methods by means of a simulation study. The simulation design is as follows: We generate data from the model $Y_{it} = m_i(\frac{t}{T}) + \beta_i X_{it} + \alpha_i + \varepsilon_{it}$, where the number of time series $i$ is set to $n = 15$ and we consider different time series lengths $T$. For simplicity, we assume that the fixed effect term $\alpha_i$ is equal to $0$ in all time series and we only include a single regressor $X_{it}$ in the model, setting the unknown parameter $\beta_i$ to $1$ for all $i$. For each $i$, the errors $\varepsilon_{it}$ follow the AR(1) model $\varepsilon_{it} = a \varepsilon_{i, t-1} + \eta_{it}$, where $a = 0.25$ and the innovations $\eta_{it}$ are i.i.d.\ normally distributed with mean $0$ and variance $0.25$. Similarly, for each $i$, the covariates $X_{it}$ follow an AR($1$) process of the form $X_{it} = a_x X_{i, t-1} + \zeta_{it}$, where $a_x = 0.5$ and the innovations $\zeta_{it}$ are i.i.d.\ normally distributed with mean $0$ and variance $1$. We assume that the covariates $X_{it}$ and the error terms $\varepsilon_{js}$ are independent from each other for all $1 \leq i,j \leq n$ and $1 \leq t, s \leq T$. To generate data under the null $H_0: m_1 = \ldots = m_n$, we let $m_i = 0$ for all $i$ without loss of generality. To produce data under the alternative, we define $m_1(u) = b \, (u - 0.5) $ with $b \in \{ 0.75, 1, 1.25 \}$ and set $m_i = 0$ for all $i \ne 1$. Hence, all trend functions are the same except for $m_1$ which is an increasing linear function. Note that the normalisation constraint $\int_0^1 m_1(u) du = 0$ is directly satisfied in this case. For each simulation exercise, we simulate $5000$ data samples.

Our multiscale test is implemented as follows: The estimators $\widehat{\beta}_i$ and $\widehat{\alpha}_i$ of the unknown parameters $\beta_i$ and $\alpha_i$ are computed as described in Section (ref). Since the errors $\varepsilon_{it}$ follow an AR($1$) process, we estimate the long-run error variance $\sigma_i^2$ by the difference-based estimator proposed in KhismatullinaVogt2020, setting the tuning parameters $q$ and $r$ to $25$ and $10$, respectively. In order to construct our test statistics, we use an Epanechnikov kernel and the grid $\mathcal{G}_T = U_T \times H_T$ with

We thus consider intervals $\mathcal{I}_{u, h} = [u-h, u+h]$ which contain $5, 15, 25, \ldots$ data points. The critical value $q_{n,T}(\alpha)$ of our test is computed by Monte Carlo methods as described in Section (ref), where we set $L=5000$.

The simulation results for our multiscale test can be found in Table (ref), which reports its actual size under $H_0$, and in Table (ref), which reports its power against different alternatives. Both the actual size and the power are computed as the number of simulations in which the test rejects the global null $H_0$ divided by the total number of simulations. Inspecting Table (ref), one can see that the actual size is fairly close to the nominal target $\alpha$ in all the considered scenarios. Hence, the test has approximately the correct size. Inspecting Table (ref), one can further see that the test has reasonable power properties. For the smallest slope $b=0.75$ and the smallest time series length $T=100$, the power is only moderate, reflecting the fact that the alternative with $b=0.75$ is not very far away from the null. However, as we increase the slope $b$ and the sample size $T$, the power increases quickly. Already for the slope value $b = 1.00$, we reach a power of $0.99$ for $T = 500$ and for all nominal sizes $\alpha$.

We next investigate the finite sample performance of the clustering algorithm from Section (ref). To do so, we consider a very simple scenario: we generate data from the model $Y_{it} = m_i(\frac{t}{T}) + \varepsilon_{it}$, that is, we assume that there are no fixed effects and no covariates. The error terms $\varepsilon_{it}$ are specified as above. Moreover, as before, we set the number of time series to $n = 15$ and we consider different time series lengths $T$. We partition the $n = 15$ time series into $N=3$ groups, each containing $5$ time series. Specifically, we set $G_1 = \{1,\ldots, 5\}$, $G_2 = \{6,\ldots, 10\}$ and $G_3 = \{11,\ldots, 15\}$, and we assume that $m_i = f_l$ for all $i \in G_l$ and all $l = 1, 2, 3$. The group-specific trend functions $f_1$, $f_2$ and $f_3$ are defined as $f_1(u) = 0$, $f_2(u) = 1 \cdot (u - 0.5)$ and $f_3(u) = (- 1) \cdot (u - 0.5)$. In order to estimate the groups $G_1$, $G_2$, $G_3$ and their number $N = 3$, we use the same implementation as before followed by the clustering procedure from Section (ref).

\addtocounter{table}{-1}

The simulation results are reported in Table (ref). The entries in Table (ref) are computed as the number of simulations for which $\widehat{N} = N$ divided by the total number of simulations. They thus specify the empirical probabilities with which the estimate $\widehat{N}$ is equal to the true number of groups $N = 3$. Analogously, the entries of Table (ref) give the empirical probabilities with which the estimated group structure $\{ \widehat{G}_1,\ldots,\widehat{G}_{\widehat{N}}\}$ equals the true one $\{G_1,G_2,G_3\}$. The results in Table (ref) nicely illustrate the theoretical properties of our clustering algorithm. According to Proposition (ref), the probability that $\widehat{N} = N$ and $\{ \widehat{G}_1,\ldots,\widehat{G}_{\widehat{N}}\} = \{G_1,G_2,G_3\}$ should be at least $(1-\alpha)$ asymptotically. For the largest sample size $T = 500$, the empirical probabilities reported in Table (ref) can indeed be seen to exceed the value $(1-\alpha)$ as predicted by Proposition (ref). For the smaller sample sizes $T=100$ and $T=250$, in contrast, the empirical probabilities are substantially smaller than $(1-\alpha)$. This reflects the asymptotic nature of Proposition (ref) and is not very surprising. It simply mirrors the fact that for the smaller sample sizes $T=100$ and $T=250$, the effective noise level in the simulated data is quite high.

Figures (ref) and (ref) give more insight into what happens for $T=100$ and $T=250$. Figure (ref) shows histograms of the $5000$ simulated values of $\widehat{N}$, while Figure (ref) depicts histograms of the number of classification errors produced by our algorithm. By the number of classification errors, we simply mean the number of incorrectly classified time series, which is formally calculated as \[ \min_{\pi \in S_{\hat{N}}} \big\{ |G_1 \setminus \widehat{G}_{\pi(1)}| +|G_2 \setminus \widehat{G}_{\pi(2)}| + |G_3 \setminus \widehat{G}_{\pi(3)}| \big\} \] with $S_{\widehat{N}}$ being the set of all permutations of $\{1, 2, \ldots, \widehat{N}\}$. The histogram of Figure (ref) for $T=100$ clearly shows that our method underestimates the number of groups ($\widehat{N} = 2$ in $4055$ cases out of $5000$). In particular, it fails to detect the difference between two out of three groups in a large number of simulations. This is reflected in the corresponding histogram of Figure (ref) which shows that there are exactly $5$ classification errors in $3924$ of the $5000$ simulation runs. In most of these cases, the estimated group structure $\{ \widehat{G}_1, \widehat{G}_{2}\}$ coincides with either $\{ G_1 \cup G_2,G_3\}$, $\{ G_1, G_2\cup G_3\}$ or $ \{ G_1 \cup G_3,G_2\}$. In summary, we can conclude that the small empirical probabilities for $T=100$ in Table (ref) are due to the algorithm underestimating the number of groups. Inspecting the histograms for $T=250$, one can see that the performance of the estimators $\widehat{N}$ and $\{ \widehat{G}_1,\ldots, \widehat{G}_{\widehat{N}} \}$ improves significantly, even though the corresponding empirical probabilities in Table (ref) are still somewhat below the target $(1-\alpha)$.

Applications

Analysis of GDP growth

In what follows, we revisit an application example from Zhang2012. The aim is to test the hypothesis of a common trend in the GDP growth time series of several OECD countries. Since we do not have access to the original dataset of Zhang2012 and we do not know the exact data specifications used there, we work with data from the following common sources: Refinitiv Datastream, the OECD.Stat database, Federal Reserve Economics Data (FRED) and the Barro-Lee Educational Attainment dataset \citep*{Barro2013}. We consider a data specification that is as close as possible to the one in Zhang2012 with one important distinction: In the original study, the authors examined $16$ OECD countries (not specifying which ones) over the time period from the fourth quarter of 1975 up to and including the third quarter of 2010, whereas we consider only $11$ countries (Australia, Austria, Canada, Finland, France, Germany, Japan, Norway, Switzerland, UK and USA) over the same time span. The reason is that we have access to data of good quality only for these 11 countries. In the following list, we specify the data for our analysis.\footnote{All data were accessed and downloaded on 7 December 2021.}

For each of the $n=11$ countries in our sample, we thus observe a quarterly time series $\mathcal{T}_i = \{(Y_{it}, \boldsymbol{X}_{it}): 1 \le t \le T \}$ of length $T = 140$, where $Y_{it} = \Delta \ln GDP_{it} := \ln GDP_{it} - \ln GDP_{i(t-1)}$ and $\boldsymbol{X}_{it} = (\Delta \ln L_{it}, \Delta \ln K_{it}, \Delta \ln H_{it})^\top$ with $\Delta \ln L_{it} := \ln L_{it} - \ln L_{i(t-1)}$, $\Delta \ln K_{it} := \ln K_{it} - \ln K_{i(t-1)}$ and $\Delta \ln H_{it} := \ln H_{it} - \ln H_{i(t-1)}$. Without loss of generality, we let $\Delta \ln GDP_{i1} = \Delta \ln L_{i1} = \Delta \ln K_{i1} = \Delta \ln H_{i1} = 0$. Each time series $\mathcal{T}_i$ is assumed to follow the model $Y_{it} = m_i(t/T) + \bm{\beta}^\top_i \boldsymbol{X}_{it} + \alpha_i + \varepsilon_{it}$, or equivalently,

for $1 \le t \le T$, where $\bm{\beta}_i = (\beta_{i, 1}, \beta_{i, 2}, \beta_{i, 3})^\top$ is a vector of unknown parameters, $m_i$ is a country-specific unknown nonparametric time trend and $\alpha_i$ is a country-specific fixed effect.

In order to test the null hypothesis $H_0: m_1 = \ldots = m_n$ with $n = 11$ in model \tagform@{(ref)}, we implement our multiscale test as follows:

Besides these choices, we construct and implement the multiscale test exactly as described in Section (ref).

The thus implemented multiscale test rejects the global null hypothesis $H_0$ at the usual significance levels $\alpha =0.01,0.05, 0.1$. This result is in line with the findings in Zhang2012 where the null hypothesis of a common trend is rejected at level $\alpha = 0.1$. The main advantage of our multiscale test over the method in Zhang2012 is that it is much more informative. In particular, it provides information about which of the $n=11$ countries have different trends and where the trends differ. This information is presented graphically in Figures (ref)--(ref). Each figure corresponds to a specific pair of countries $(i, j)$ and is divided into three panels (a)--(c). Panel (a) shows the augmented time series $\{\widehat{Y}_{it}: 1 \le t \le T\}$ and $\{\widehat{Y}_{jt}: 1 \le t \le T\}$ for the two countries $i$ and $j$ that are compared. Panel (b) presents smoothed versions of the time series from (a), in particular, it shows local linear kernel estimates of the two trend functions $m_i$ and $m_j$, where the bandwidth is set to $14$ quarters (that is, to $0.1$ in terms of rescaled time) and an Epanechnikov kernel is used. Panel (c) presents the results produced by our test for the significance level $\alpha = 0.05$: it depicts in grey the set $\mathcal{S}^{[i, j]}(\alpha)$ of all the intervals for which the test rejects the local null $H_0^{[i, j]}(u, h)$. The set of minimal intervals $\mathcal{S}^{[i, j]}_{min}(\alpha) \subseteq \mathcal{S}^{[i, j]}(\alpha)$ is highlighted in black. According to \tagform@{(ref)}, we can make the following simultaneous confidence statement about the intervals plotted in panels (c) of Figures (ref)--(ref): we can claim, with confidence of about $95\%$, that there is a difference between the functions $m_i$ and $m_j$ on each of these intervals.

Out of $55$ pairwise comparisons, our test detects differences for $11$ pairs of countries $(i,j)$. These $11$ cases are presented in Figures (ref)--(ref). In $9$ cases (Figures (ref)--(ref)), one of the involved countries is Norway. Inspecting the trend estimates in panels (b) of Figures (ref)--(ref), the Norwegian trend estimate can be seen to exhibit a strong downward movement at the end of the observation period, whereas the other trend estimates show a much less pronounced downward movement (or even a slight upward movement). According to our test, this is a significant difference between the Norwegian and the other trend functions rather than an artefact of the sampling noise: In all $9$ cases, the test rejects the local null for at least one interval which covers the last $10$ years of the analysed time period (from the first quarter in $2000$ up to the third quarter in $2010$). Apart from these differences at the end of the sampling period, our test also finds differences in the beginning, however, only for part of the pairwise comparisons.

Figures (ref) and (ref) present the results of the pairwise comparison between Australia and France and between the USA and France, respectively. In both cases, our test detects differences between the GDP trends only in the beginning of the considered time period. In the case of Australia and France, it is clearly visible in the raw data (panel (a) in Figure (ref)) that there is a difference between the trends, whereas this is not so obvious in the case of the USA and France. According to our test, there are indeed significant differences in both cases. In particular, we can claim with confidence at least 95%, that there are differences between the trends of the USA and France (of Australia and France) up to the fourth quarter in 1991 (the fourth quarter in 1986), but there is no evidence of any differences between the trends from 1992 (1987) onwards.

We next apply our clustering techniques to find groups of countries that have the same time trend. We implement our HAC algorithm with $\alpha = 0.05$ and the same choices as detailed above. The dendrogram that depicts the clustering results is plotted in Figure (ref). The number of clusters is estimated to be $\widehat{N} = 3$. The rectangles in Figure (ref) indicate the $\widehat{N} = 3$ clusters. In particular, each rectangle is drawn around the branches of the dendrogram that correspond to one of the three clusters. Figure (ref) depicts local linear kernel estimates of the $n=11$ GDP time trends (calculated from the augmented time series $\widehat{Y}_{it}$ with bandwidth $0.1$ and Epanechnikov kernel). Their colour indicates which cluster they belong to.

The results in Figures (ref) and (ref) show that there is one cluster which consists only of Norway (plotted in red). As we have already discussed above and as becomes apparent from Figure (ref), the Norwegian trend exhibits a strong downward movement at the end of the sampling period, whereas the other trends show a much more moderate downward movement (if at all). This is presumably the reason why the clustering procedure puts Norway in a separate cluster. The algorithm further finds two other clusters, one consisting of the 5 countries Australia, Finland, Germany, Japan and the USA (plotted in blue in Figures (ref) and (ref)) and the other one consisting of the 5 countries Austria, Canada, France, Switzerland and the UK (plotted in green in Figures (ref) and (ref)). Visual inspection of the trend estimates in Figure (ref) suggests that the GDP time trends in the blue cluster exhibit more pronounced decreases and increases than the GDP time trends in the green cluster. Hence, overall, the clustering procedure appears to produce a reasonable grouping of the GDP trends.

Analysis of house prices

We next analyse a historical dataset on nominal annual house prices from Knoll2017 that contains data for $14$ advanced economies covering $143$ years from $1870$ to $2012$. In our analysis, we consider 8 countries (Australia, Belgium, Denmark, France, Netherlands, Norway, Sweden and USA) over the time period 1890--2012. The data for all these countries except one (Belgium) contain no missing values, and for Belgium there are only five missing observations\footnote{The missing values in the Belgium time series span five years during World War I.} which we impute by linear interpolation. The time series of the other $6$ countries contain more than $10$ missing values each, which is why we exclude them from our analysis.

We deflate the nominal house prices with the corresponding consumer price index (CPI) to obtain real house prices ($HP$). Variables that can potentially influence the average house prices are numerous, and there seems to be no general consensus about what the main determinants are. Possible determinants include, but are not limited to, demographic factors such as population growth (Holly2010, Wang2014, Churchill2021); fundamental economic factors such as real GDP (Huang2013, Churchill2021), interest rate and inflation (Abelson2005, Otto2007, Huang2013, Jorda2015); urbanisation (Chen2011, Wang2017); government subsidies and regulations (Malpezzi1999); stock markets (Gallin2006); etc. In our analysis, we focus on the following determinants of the average house prices: real GDP ($GDP$), population size ($POP$), long-term interest rate ($I$) and inflation ($INFL$) which is measured as change in CPI. Most other factors (such as government regulations, construction costs, and real wages) vary rather slowly over time and can be captured by time trend, fixed effects and slope heterogeneity. Data for CPI, real GDP, population size and long-term interest rate are taken from the Jordà-Schularick-Taylor Macrohistory Database\footnote{See Jorda2017 for a detailed description of the variable construction.}, which is freely available at \url{http://www.macrohistory.net/data/} (accessed on 13 January 2022).

In summary, we observe a panel of $n = 8$ time series $\mathcal{T}_i = \{(Y_{it}, \boldsymbol{X}_{it}): 1 \le t \le T \}$ of length $T = 123$ for each country $i \in \{1,\ldots, 8\}$, where $Y_{it} = \ln HP_{it}$ and $\boldsymbol{X}_{it} = (\ln GDP_{it}, \ln POP_{it}, I_{it}, INFL_{it})^\top$. For each $i$, the time series $\mathcal{T}_i$ is assumed to follow the model $Y_{it} = m_i(t/T) + \bm{\beta}^\top_i \boldsymbol{X}_{it} + \alpha_i + \varepsilon_{it}$, or equivalently,

for $1 \le t \le T$, where $\bm{\beta}_i = (\beta_{i, 1}, \beta_{i, 2}, \beta_{i, 3}, \beta_{i, 4})^\top$ is a vector of unknown parameters, $m_i$ is a country-specific unknown nonparametric time trend and $\alpha_i$ is a fixed effect.

The inclusion of a nonparametric trend function $m_i$ in model \tagform@{(ref)} is supported by the literature. Ugarte2009, for example, model the trend in average Spanish house prices by means of splines. Winter2022 include a long-term stochastic trend component when describing the dynamic behaviour of real house prices in $8$ advanced economies. Including a nonparametric trend function when modelling the evolution of house prices is also the main conclusion in Zhang2016, where it is shown that the time series of logarithmic US house prices is trend-stationary, i.e., can be transformed into a stationary series by subtracting a deterministic trend.

We implement the multiscale test from Section (ref) in the same way as in the previous application example with one minor modification: we let $\mathcal{G}_T = U_T \times H_T$ with

We thus take into account all locations $u$ on an equidistant grid $U_T$ with step length $1/T$ and all scales $h=2/T, 7/T, 12/T,\ldots$ with $\log T /T \le h \le 1/4$. This implies that each interval $\mathcal{I}_{u, h} =[u-h,u+h]$ with $(u,h) \in \mathcal{G}_T$ spans $5, 15, 25, \ldots$ years. The lower bound $\log T / T$ is motivated by Assumption (ref). As in Section (ref), we assume that for each $i$, the error process $\mathcal{E}_i = \{\varepsilon_{it}: 1 \leq t \leq T\}$ follows an AR($p_i$) model and we estimate the long-run variances $\sigma_i^2$ by the difference-based estimator from KhismatullinaVogt2020 with tuning parameters $q$ and $r$ equal to $15$ and $10$, respectively. We choose $p_i$ by minimizing the BIC. For $7$ out of $8$ countries the order $p_i$ determined by BIC\footnote{Applying other information criteria such as FPE, AIC and HQ yields exactly the same results in these cases.} is equal to $1$. For the sake of simplicity, we thus assume that $p_i = 1$ for all $i$.

We are now ready to apply our test to the data. The overall null hypothesis $H_0$ is rejected at levels $\alpha = 0.05$ and $\alpha = 0.1$. The detailed test results for the significance level $\alpha =0.05$ are presented in Figures (ref)--(ref). As in Section (ref), each figure corresponds to the comparison of a pair of countries $(i, j)$ for which our test detects differences between the trends. Panel (a) shows the augmented time series $\{\widehat{Y}_{it}: 1 \le t \le T\}$ and $\{\widehat{Y}_{jt}: 1 \le t \le T\}$ for the two countries $i$ and $j$ under consideration. Panel (b) presents smoothed versions of the time series from (a) with the bandwidth window covering $15$ years. Panel (c) presents the test results for the level $\alpha=0.05$. As before, the set of intervals $\mathcal{S}^{[i, j]}(\alpha)$ for which our test rejects and the set of minimal intervals $\mathcal{S}^{[i, j]}_{min}(\alpha) \subseteq \mathcal{S}^{[i, j]}(\alpha)$ are depicted in grey and black, respectively. According to \tagform@{(ref)}, we can make the following simultaneous confidence statement about the intervals plotted in panels (c) of Figures (ref)--(ref): we can claim, with confidence of about $95\%$, that there is a difference between the trends $m_i$ and $m_j$ on each of these intervals.

Overall, our findings are in line with the observations in Knoll2017, where the authors find considerable cross-country heterogeneity in the house price trends. The authors note that before World War II, the countries exhibit similar trends in real house prices, while the trends start to diverge sometime after World War II. This fits with our findings in Figures (ref) and (ref) which show that our test detects differences between the trends of Australia and the Netherlands (of Belgium and the Netherlands) starting only from 1968 (1966) onwards. Contrary to Knoll2017, however, our test also finds significant differences in the first half of the observed time period, specifically, between the time trends of Belgium and Denmark and of Belgium and the USA (Figures (ref) and (ref), respectively). This discrepancy in the results is most certainly due to the fact that the method used in Knoll2017 does not account for effects of other factors such as GDP or population growth, while our test allows us to include various determinants of the average house prices in model \tagform@{(ref)}.

We next apply the clustering procedure from Section (ref) to the data. As in the previous application example, we set $\alpha = 0.05$. The results are displayed in Figures (ref) and (ref). Specifically, Figure (ref) shows the dendrogram with the results of the HAC algorithm and Figure (ref) presents local linear kernel estimates of the trend curves (calculated with a bandwidth window of $15$ years and an Epanechnikov kernel). The number of clusters is estimated to be $\widehat{N} = 3$. As before, the coloured rectangles in Figure (ref) are drawn around the countries that belong to the same cluster and the same colours are used to display the trend estimates in Figure (ref).

Inspecting the results, we can see that there is one cluster consisting only of Belgium (plotted in green). Figure (ref) suggests that the Belgium time trend indeed evolves somewhat differently from the other trends in the first $30$ years of the observed time period. The algorithm further detects a cluster that consists of two countries: France and the Netherlands (plotted in blue). The time trends of these two countries display some kind of dip around 1950, which is not present in the time trends of the other countries. Overall, our algorithm thus appears to produce a reasonable clustering of the house price trends.

Acknowledgements

Financial support by the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation), Germany – grant VO 2503/1-1, project number 430668955 – is gratefully acknowledged.

{ {0.55em} }

\allowdisplaybreaks[3]