Distributed additive hazards regression analysis of multi-site current status data without using individual-level data

In multi-site studies, sharing individual-level information across multiple data-contributing sites usually poses a significant risk to data security. Thus, due to privacy constraints, analytical tools without using individual-level data have drawn considerable attention to researchers in recent years. In this work, we consider regression analysis of current status data arising from multi-site cross-sectional studies and develop two distributed estimation methods tailored to the unstratified and stratified additive hazards models, respectively. In particular, instead of utilizing the individual-level data, the proposed methods only require transferring the summary statistics from each site to the analysis center, which achieves the aim of privacy protection. We establish the asymptotic properties of the proposed estimators, including the consistency and asymptotic normality. Specifically, the distributed estimators derived are shown to be asymptotically equivalent to those based on the pooled individual-level data. Simulation studies and an application to a multi-site gonorrhea infection data set demonstrate the proposed methods’ satisfactory performance and practical utility.

The data that support the findings in this paper are available on request from the corresponding author. The data are not publicly available due to privacy.

We are grateful to the associate editor and two reviewers for their insightful comments and suggestions that greatly improve this article. Shuwei Li’ research was partially supported by the Nature Science Foundation of Guangdong Province of China (Grant No. 2022A1515011901) and National Statistical Science Research Project (Grant No. 2022LY041). Xinyuan Song’ research was partially supported by GRF Grant (Grant No. 14303622) from the Research Grant Council of the Hong Kong Special Administrative Region.

Authors and Affiliations

1. School of Economics and Statistics, Guangzhou University, Guangzhou, China

   Peiyao Huang & Shuwei Li

2. Department of Statistics, Chinese University of Hong Kong, Sha Tin, Hong Kong, China

   Xinyuan Song

Contributions

Conceptualization, SL and XS; methodology, PH, SL and XS; software, PH; resources, SL; data curation, SL; writing—original draft preparation, PH, SL and XS; supervision, XS; funding acquisition, SL and XS; All authors reviewed the manuscript.

Corresponding authors

Correspondence to Shuwei Li or Xinyuan Song.

Conflict of interest

The authors declare that they have no conflict of interest.

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Appendix A: Regularity conditions

Recall that \(\Omega = \{1, \dots , n\}\) is the index set of the observed data from K data-contributing sites. \(\Omega _k = \{i: i\) is from site k, for \(i = 1, \ldots ,n\}\) denotes the index set of the kth site, and \(n_k\) is the size of \(\Omega _k\), then \(n=n_1+ \ldots +n_K\). The size of \(\widetilde{\Omega }\) is d, where \(\widetilde{\Omega }=\{i:\delta _i =1,\) for \(i=1, \ldots , n\}\). For each \(k = 1, \ldots , K\), \(\widetilde{\Omega }_k=\{i: \delta _i=1, i\) is from site k,  for \(i=1, \ldots , n\}\). Denoted by \(\varvec{\beta }_0\) the true value of \(\varvec{\beta }\) related to the population formed by the K sites. Denoted by \(\varvec{\beta }_{0(k)}\) the true value of \(\varvec{\beta }\) related to the kth site. Let \(\varvec{\theta _0}=(\varvec{\beta }_{0(1)}^{\top },\ldots ,\varvec{\beta }_{0(K)}^{\top },\varvec{\beta }_{0}^{\top })^{\top }\) denote the true value of \(\varvec{\theta }=(\varvec{\beta }_{(1)}^{\top },\ldots ,\varvec{\beta }_{(K)}^{\top },\varvec{\beta }^{ \top })^{\top }\).

Let \(S^{(m)}(\varvec{\beta },t)=\sum _{j \in \Omega } Y_{j}(t)e^{-\varvec{\beta }^{\top } \varvec{Z}_{j}t}{(\varvec{Z}_{j}t)}^{\otimes m} \), \(S_k^{(m)}(\varvec{\beta },t)=\sum _{j \in \Omega _k} Y_{j(k)}(t)e^{-\varvec{\beta }_{(k)}^{\top } \varvec{Z}_{j(k)}t}{(\varvec{Z}_{j(k)}t)}^{\otimes m}, \) for \(k=1, \ldots , K\) and \(m = 0,1,2\), where \(\varvec{a}^{\otimes 0} = 1, \varvec{a}^{\otimes 1} = \varvec{a}\), and \( \varvec{a}^{\otimes 2} = \varvec{a}\varvec{a}^\top \) for a column vector \(\varvec{a}\). Define \(E_k(\varvec{\beta },t)= S_k^{(1)}(\varvec{\beta },t)/S_k^{(0)}(\varvec{\beta },t)\), \(V_k(\varvec{\beta },t)= S_k^{(2)}(\varvec{\beta },t)/S_k^{(0)}(\varvec{\beta },t) - E_k(\varvec{\beta },t)^{\otimes {2}}\) for \(t \in (0, \tau ]\) and \(\varvec{\beta }\in \mathcal {B}\), where \(\mathcal {B}\) is a compact set in \(\mathbb {R}^p\). For a matrix \(\varvec{A}\) or a vector \(\varvec{a}\), \(||\varvec{A}||=\sup _{i,j}|a_{ij}|\) and \(||\varvec{a}||=\sup _{i}|a_{i}|\), where \(a_{ij}\) is the (i, j)th element of \(\varvec{A}\) and \(a_{i}\) is the ith component of \(\varvec{a}\). Let “\(\overset{P}{\rightarrow }\)” and “\(\overset{D}{\rightarrow }\)” denote the convergence in probability and the convergence in distribution, respectively.

To establish the asymptotic properties of proposed estimator \(\hat{\varvec{\beta }}\), we need the following regularity conditions:

1. (C1)

   \(P(Y(t) = 1 \mid \varvec{Z}) > 0\) for all \(t \in (0, \tau ]\), where \(\tau \) satisfies \(P(C \ge \tau ) > 0\) and \(Y(t) = I(C \ge t)\).

2. (C2)

   \(\varvec{Z}_i\) is bounded for \(i=1,2,\ldots ,n\). The true value of H(t) denoted by \(H_0(t)\) at \(\tau \) is finite. For \(k = 1, \ldots , K\), the true value of \(H_k(t)\) denoted by \(H_{0(k)}(t)\) at \(\tau \) is finite, where \(\textrm{d} H_k(t) = e^{-\Lambda _k (t)}\textrm{d}\Lambda _c(t)\) and \(\Lambda _k (t)=\int _{0}^{t}\lambda _k (s)\textrm{d}s\).

3. (C3)

   The true regression vector \(\varvec{\beta }_0\) lies in the interior of \(\mathcal {B}\) and there exist functions \(s^{(m)}(\varvec{\beta },t)\) with \(m=0,1,2\) defined on \(\mathcal {B}\times (0,\tau ]\) satisfying

   1. (a)

      \(\sup \limits _{\varvec{\beta }\in \mathcal {B}, t\in (0,\tau ]} ||n^{-1}S^{(m)}(\varvec{\beta },t)-s^{(m)}(\varvec{\beta },t)|| \overset{P}{\rightarrow }\textbf{0}\),   as \(n\rightarrow \infty \).

   2. (b)

      \(s^{(0)}(\varvec{\beta },t)\) is bounded away from 0.

   3. (c)

      For \(m=0,1,2\), \(s^{(m)}(\varvec{\beta },t)\) is a uniformly continuous function of \(\varvec{\beta }\) in \((0, \tau ]\), where \(s^{(1)}(\varvec{\beta },t) = \partial {s^{(0)}(\varvec{\beta },t)}/{\partial \varvec{\beta }} \) and \(s^{(2)}(\varvec{\beta },t)=\partial ^2{s^{(0)}(\varvec{\beta },t)}/{\partial \varvec{\beta }\partial \varvec{\beta }^{\top }}.\)

   4. (d)

      For \(\varvec{\beta }\in \mathcal {B}\), let \(e(\varvec{\beta },t)=s^{(1)}(\varvec{\beta },t)/s^{(0)}(\varvec{\beta },t)\), \(v(\varvec{\beta },t)=s^{(2)}(\varvec{\beta },t)/s^{(0)}(\varvec{\beta },t)-e(\varvec{\beta },t) e(\varvec{\beta },t)^{\top }.\) \(\varvec{\Sigma }^*(\varvec{\beta }_0)=\int _0^{\tau }v(\varvec{\beta }_0,u)s^{(0)}(\varvec{\beta }_0,u)\textrm{d} H_0(u)\) is positive definite.

4. (C4)
   1. (a)

      Underlying parameters are identical across all sites, that is, \(\varvec{\beta }_{0(k)}=\varvec{\beta }_0\) for \(k=1,\ldots , K\);

   2. (b)

      Limiting processes of \(S_{k}^{(m)}(\varvec{\beta },t)\)’s are identical across all sites, that is, for \(k=1,\ldots , K\), \(\sup \limits _{\varvec{\beta }\in \mathcal {B}, t\in (0,\tau ]} ||n_k^{-1}S_k^{(m)}(\varvec{\beta },t)-s^{(m)} (\varvec{\beta },t)|| \overset{P}{\rightarrow }\ 0\), as \(n_k\rightarrow \infty \) for \(\varvec{\beta }\in \mathcal {B}\) and \(t \in (0, \tau ]\).

Conditions (C1)–(C3) are standard assumptions in the literature of survival analysis, see for instance Andersen and Gill (1982). Condition (C4) is usually referred to as the homogeneity assumption, which implies that data from different sites follow the same underlying model or different models but with same regression parameters (Li et al. 2023, for example).

Appendix B: Proof of Theorem 1

Let \(\hat{\varvec{\beta }}\) denote the estimator of \(\varvec{\beta }\) by solving \(\widetilde{\varvec{U}}(\varvec{\beta })=0\). For each \(k = 1, \dots , K\), we first obtain the estimate of \(\varvec{\beta }\) denoted by \(\hat{\varvec{\beta }}_{(k)}\) by solving \(\varvec{U}_k^*(\varvec{\beta }) = 0\) within site k, where \(\varvec{U}^*_k(\varvec{\beta })\) has the same form as \(\varvec{U}^* (\varvec{\beta })\) (3) after replacing \(\Omega \) with \(\Omega _k\). Define \(\hat{\varvec{\theta }} = (\hat{\varvec{\beta }}_{(1)}^{\top }, \ldots , \hat{\varvec{\beta }}_{(K)}^{\top }, \hat{\varvec{\beta }}^{\top })^{\top }\). By applying Taylor series expansion about \(\varvec{\beta }_{0}\) to \(\varvec{U}_k^*(\hat{\varvec{\beta }}_{(k)}) \), we have

where \(\varvec{\beta }^*_{(k)}\) is between \(\hat{\varvec{\beta }}_{(k)}\) and \(\varvec{\beta }_0\) and \(\varvec{I}^*_k(\varvec{\beta })\) has the same form as \(\varvec{I}^*(\varvec{\beta })\) (4) after replacing \(\Omega \) with \(\Omega _k\). By the arguments given in Andersen and Gill (1982) and Conditions (C1)–(C4), it can be shown that \(\hat{\varvec{\beta }}_{(k)} \overset{P}{\rightarrow }\ \varvec{\beta }_0\), and \(n_k^{-1}\varvec{I}_k^*(\varvec{\beta }^*_{(k)})\overset{P}{\rightarrow }\varvec{\Sigma }^*(\varvec{\beta }_0)\) for any \(\varvec{\beta }_{(k)}^*\) that converges in probability to \(\varvec{\beta }_0\). For each \(k=1, \ldots , K\), we follow Lin et al. (1998) and can conclude that \(n_k^{-1/2} \varvec{U}_k^*(\varvec{\beta }_0) \overset{D}{\rightarrow } \varvec{N}(\varvec{0},\varvec{\Sigma }^*(\varvec{\beta }_0))\) and \(n_k^{1/2}(\hat{\varvec{\beta }}_{(k)}-\varvec{\beta }_0)\overset{D}{\rightarrow } N(\varvec{0}, \varvec{\Sigma }^*(\varvec{\beta }_0)^{-1})\) under Conditions (C1)–(C4). Note that \(\varvec{\Sigma }^*(\varvec{\beta }_0) = \lim \limits _{n_k \rightarrow \infty }n_k^{-1}\)\(\varvec{I}_k^*(\varvec{\beta }_0)\), and under Condition (C4), we also have \(\varvec{\Sigma }^*(\varvec{\beta }_0) = \lim \limits _{n \rightarrow \infty }n^{-1}\varvec{I}(\varvec{\beta }_0)\).

Consider the following Taylor series expansion

Then \(S^{(0)}(\varvec{\beta },t)\), \(S^{(1)}(\varvec{\beta },t)\) and \(S^{(2)}(\varvec{\beta },t)\) can be approximated by

and

respectively. Under Conditions (C1)–(C4), we have \(\sup \limits _{\varvec{\beta }\in \mathcal {B}, t\in (0,\tau ]}\)\( ||n^{-1}\widetilde{S}^{(m)}(\varvec{\beta },t)-s^{(m)} (\varvec{\beta },t)|| \overset{P}{\rightarrow }\ 0\) as \(n\rightarrow \infty \), for \(m=0,1,2\) and \(\hat{\varvec{\beta }}_{(k)}\) that is close to \(\varvec{\beta }\). The pseudo log-likelihood function related to the score function (7) is given by

The first and second derivatives of \(\widetilde{l}(\varvec{\beta })\) with respect to \(\varvec{\beta }\) are

and

respectively. By following Andersen and Gill (1982) and Kalbfleisch and Prentice (2002), we define

for \(t \le \tau \), where \(\widetilde{l}(\varvec{\beta }, t)=-\sum _{i \in \widetilde{\Omega }} \int _0^{t} \{\varvec{\beta }^{\top }\varvec{Z}_{i}t + \textrm{log}\)\( [\widetilde{S}^{(0)}(\varvec{\beta },t) ] \} \textrm{d}N_{i}(u)\). Note that \(\hat{\varvec{\beta }}\) maximizes \(X(\varvec{\beta },\tau )\), we can obtain the compensator of \(X(\varvec{\beta },t)\), which is given by

By following Andersen and Gill (1982), it can be shown that \(X(\varvec{\beta }, t)-\tilde{X}(\varvec{\beta }, t)\) converges to 0 for any \(t \in (0,\tau ]\). Define

Under Conditions (C1)–(C4), it can be concluded that \(\tilde{X}(\varvec{\beta },\tau )\overset{P}{\rightarrow }f(\varvec{\beta })\). By applying the Lenglart’s inequality given in Andersen and Gill (1982), we have \(X(\varvec{\beta },\tau )\overset{P}{\rightarrow }f(\varvec{\beta })\). The first and second derivatives of \(f(\varvec{\beta })\) with respect to \(\varvec{\beta }\) are

and

respectively. Note that \(\partial f(\varvec{\beta })/\partial \varvec{\beta }|_{\varvec{\beta }=\varvec{\beta }_0} = 0\) and the negative of \(\partial ^{2} f(\varvec{\beta })/\partial \varvec{\beta }\partial \varvec{\beta }^{\top }\) is a positive definite matrix, \(\varvec{\beta }_0\) is the unique maximum of \(f(\varvec{\beta })\). Therefore, by following Andersen and Gill (1982), we can conclude that \(\hat{\varvec{\beta }} \overset{P}{\rightarrow }\ \varvec{\beta }_0\) as \(n {\rightarrow } \infty \).

Notice that

By applying Taylor series expansion about \(\varvec{\theta }_0\) to \(\widetilde{\varvec{U}}(\hat{\varvec{\theta }})\) and following the arguments given in Li et al. (2023), we have

where \(\varvec{\theta ^*}\) is between \(\hat{\varvec{\theta }}\) and \(\varvec{\theta }_0\), \(\widetilde{\varvec{I}}_{k} (\varvec{\theta })= \partial \widetilde{\varvec{U}}(\varvec{\theta })/\partial \varvec{\beta }_{(k)}^{\top }\) and \(\widetilde{\varvec{I}}_{K+1} (\varvec{\theta })= \partial \widetilde{\varvec{U}}(\varvec{\theta })/\partial \varvec{\beta }^{\top }\). The consistency of \(\hat{\varvec{\theta }}\) indicates that \(n_k^{-1}\widetilde{\varvec{I}}_{k}({\varvec{\theta }^*}) \overset{P}{\rightarrow } \varvec{\Sigma }_{k}(\varvec{\theta }_0)\) and \(n^{-1}\widetilde{\varvec{I}}_{K+1}({\varvec{\theta }^*}) \overset{P}{\rightarrow } \varvec{\Sigma }_{K+1}(\varvec{\theta }_0)\) for any \(\varvec{\theta }^*\) that converges in probability to \(\varvec{\theta }_0\), where \(\varvec{\Sigma }_{k}(\varvec{\theta }_0) = \lim \limits _{n_k \rightarrow \infty }n_k^{-1}\widetilde{\varvec{I}}_k(\varvec{\theta }_0)\) and \(\varvec{\Sigma }_{K+1}(\varvec{\theta }_0) = \lim \limits _{n \rightarrow \infty }n^{-1}\widetilde{\varvec{I}}_{K+1}(\varvec{\theta }_0)\).

Next, we show that \(\textrm{Cov}(\varvec{U}^*_k({\varvec{\beta }_{0}}),\widetilde{\varvec{U}}({\varvec{\theta }_0}))=\varvec{0}\).

For \(k=1,\ldots ,K\), the observed data from the kth site are denoted by \(\mathcal {O} = \{ (C_j, \delta _j, \varvec{Z}_j), j\in \Omega _k \}\). Without loss of generality, we temporarily assume \(K=2\), the covariance between \(\varvec{U}_1^{*}(\varvec{\beta }_{0})\) and \(\widetilde{\varvec{U}}(\varvec{\theta }_0)\) is

Since \(\mathbb {E}(\varvec{U}_1^*(\varvec{\beta }_{0}))=\varvec{0}\), it can be conclude that \(\textrm{Cov}(\varvec{U}_1^*(\varvec{\beta }_{0}),\)\( \widetilde{\varvec{U}}(\varvec{\theta }_0))=\varvec{0}\). Similarly, we have \(\textrm{Cov}(\varvec{U}_2^*(\varvec{\beta }_{0}),\widetilde{\varvec{U}}(\varvec{\theta }_0))=\varvec{0}\). Additionally, owing to the independence of K data-contributing sites, we have \( \textrm{Cov}(\varvec{U}_i^*(\varvec{\beta }_{0}), \varvec{U}_j^*(\varvec{\beta }_{0}))=\varvec{0}\) for \(i \ne j\). By applying the Central Limit Thoerem (Andersen and Gill 1982; Kalbfleisch and Prentice 2002), we can conclude that

where

Furthermore, by (13) and (14), we have

where

Hence, we have

as \(n {\rightarrow } \infty \). Then by (15), we can conclude that \(n^{1/2}(\hat{\varvec{\beta }}-\varvec{\beta }_0)\) converges in distribution to a zero-mean normal random vector with covariance matrix

which corresponds to the bottom right corner block of \( \textbf{P}^{-1}\textbf{H}(\textbf{P}^{\top })^{-1}\). Note that

where

and

where

and

In (14), under Conditions (C1)–(C4), we know that \(\varvec{\Sigma }_{K+1}(\varvec{\theta }_0)\) equals the limit of \(n^{-1}\widetilde{\varvec{I}}_{K+1}(\varvec{\theta }_0)\), \(\varvec{\Sigma }_k(\varvec{\theta }_0)\) equals the limit of \(n_k^{-1}\widetilde{\varvec{I}}_{k}(\varvec{\theta }_0)\) and \(\varvec{\Sigma }^*(\varvec{\beta }_{0})\) is the limit of \(n^{-1}{\varvec{I}} (\varvec{\beta }_0)\). Therefore, since \(\varvec{\Sigma }_k(\varvec{\theta }_0) =\varvec{0}\) for \(k=1,\ldots ,K\), \(\varvec{\Sigma }_{K+1}(\varvec{\theta }_0)=\varvec{\Sigma }^*(\varvec{\beta }_{0})\), we can conclude from (16) that

Notice that when the pooled individual-level data from the K sites are accessible, the estimate of \(\varvec{\beta }\) can be obtained by solving \(\varvec{U}(\varvec{\beta })=0\). Based on the arguments provided in Andersen and Gill (1982) and Conditions (C1)–(C4), the limiting covariance matrix of the estimator of \(\varvec{\beta }\) based on the pooled individual-level data is the same as the proposed distributed estimator \(\hat{\varvec{\beta }}\). This implies that the distributed estimator \(\hat{\varvec{\beta }}\) is asymptotically equivalent to the estimator obtained from analyzing the pooled individual-level data.

By the conclusions obtained above, we know that \(\varvec{\Sigma }_{K+1}(\varvec{\theta _0})\), \(\varvec{\Sigma }_{k}(\varvec{\theta _0})\) and \(\varvec{\Sigma }^{*}(\varvec{\theta _0})\) can be consistently estimated by \(n^{-1}\widetilde{\varvec{I}}( \hat{\varvec{\theta }})\), \(n_k^{-1}\varvec{\widetilde{I}}_{k}(\hat{\varvec{\theta }})\) and \(n_k^{-1}\varvec{I}_k^*(\hat{\varvec{\beta }}_{(k)})\), respectively. In light of (16), we can obtain a consistent estimate of the covariance matrix of \(\hat{\varvec{\beta }}\), which is given by

Appendix C: Proof of Theorem 2

Under Conditions (C1)–(C4) and by following Andersen and Gill (1982), we have \(\hat{\varvec{\beta }}_{(k)}\overset{P}{\rightarrow }\ \varvec{\beta }_0\), \(n_k^{-1}\varvec{I}^*_k(\hat{\varvec{\beta }}_{(k)})\overset{P}{\rightarrow } \varvec{\Sigma }^*(\varvec{\beta }_0)\), \(n_k^{-1}\varvec{\Phi }_k(\hat{\varvec{\beta }}_{(k)})\overset{D}{\rightarrow } \varvec{N}(\varvec{0},\varvec{\Sigma }^*(\varvec{\beta }_0))\) and \(n_k^{1/2}(\hat{\varvec{\beta }}_{(k)}-\varvec{\beta }_0)\overset{D}{\rightarrow } N(\varvec{0}, \varvec{\Sigma }^*(\varvec{\beta }_0)^{-1})\) as \(n_k {\rightarrow } \infty \), where \(\varvec{\Sigma }^*(\varvec{\beta }_0) = \lim \limits _{n_k \rightarrow \infty }n_k^{-1}\varvec{I}_k^*(\varvec{\beta }_0)\). Let \(\hat{\varvec{\beta }}\) denote the estimator of \(\varvec{\beta }\) by solving \(\widetilde{\varvec{\Phi }}(\varvec{\beta })=0\). Under Condition (C4), we also have that \(\varvec{\Sigma }^*(\varvec{\beta }_0) = \lim \limits _{n \rightarrow \infty }n^{-1}\varvec{H}(\varvec{\beta }_0)\). Notice that \(S_k^{(0)}(\varvec{\beta })\), \(S_k^{(1)}(\varvec{\beta })\) and \(S_k^{(2)}(\varvec{\beta })\) are approximated by

and

respectively. The pseudo log-likelihood function related to the score function (11) is given by

The first and second derivatives of \(\widetilde{l}(\varvec{\beta })\) with respect to \(\varvec{\beta }\) are

and

respectively. By following the similar arguments given in Appendix B, we can conclude that \(\hat{\varvec{\beta }}\overset{P}{\rightarrow } \varvec{\beta }_0\) as \(n {\rightarrow } \infty \). Notice that

To derive the asymptotic normality of \(\hat{\varvec{\beta }}\), we first apply Taylor series expansion about \(\varvec{\theta }_0\) to \(\widetilde{\varvec{\Phi }}(\hat{\varvec{\theta }})\) and obtain

where \(\varvec{\theta ^*}\) is between \(\hat{\varvec{\theta }}\) and \(\varvec{\theta }_0\), \(\varvec{\widetilde{H}}_{k} (\varvec{\theta })= \partial \varvec{\widetilde{\Phi }}(\varvec{\theta })/\partial \varvec{\beta }_{(k)}^{\top }\) and \(\varvec{ \widetilde{H}}_{K+1} (\varvec{\theta })= \partial \varvec{\widetilde{\Phi }}(\varvec{\theta })/\partial \varvec{\beta }^{\top }\). By following Lin et al. (1998) and Li et al. (2023), for \(k=1,\ldots K\), it can be concluded that \(n_k^{-1}\varvec{\widetilde{H}}_{k}({\varvec{\theta }^*}) \overset{P}{\rightarrow } \varvec{\Sigma }_{k}(\varvec{\theta }_0)\) and \(n^{-1}\varvec{ \widetilde{H}}_{K+1}({\varvec{\theta }^*}) \overset{P}{\rightarrow } \varvec{\Sigma }_{K+1}(\varvec{\theta }_0)\) for any \(\varvec{\theta }^*\) that converges in probability to \(\varvec{\theta }_0\), where \(\varvec{\Sigma }_{k}(\varvec{\theta }_0) = \lim \limits _{n_k \rightarrow \infty }n_k^{-1}\varvec{\widetilde{H}}_k(\varvec{\theta }_0)\) and \(\varvec{\Sigma }_{K+1}(\varvec{\theta }_0) = \lim \limits _{n \rightarrow \infty }n^{-1}\varvec{\widetilde{H}}_{K+1}(\varvec{\theta }_0)\). By following the similar arguments given in Appendix B, it can be concluded that

where

Furthermore, by Taylor expansion about \(\varvec{\beta }_0\) to \(\varvec{\Phi }_k(\hat{\varvec{\beta }}_{(k)})\) and the equation (17), we have

where

Therefore, we have that

as \(n {\rightarrow } \infty \). Then we can conclude that \(n^{1/2}(\hat{\varvec{\beta }}-\varvec{\beta }_0)\) converges in distribution to a zero-mean normal random vector with covariance matrix

which corresponds to the bottom right corner block of \( \textbf{P}^{-1}\textbf{H}(\textbf{P}^{\top })^{-1}\). Note that

where

and

where

and

In (17), by Conditions (C1)–(C4), we know that \(\varvec{\Sigma }_{K+1}(\varvec{\theta }_0)\) is the limit of \(n^{-1}\varvec{\widetilde{H}}_{K+1}(\varvec{\theta }_0)\), \(\varvec{\Sigma }_k(\varvec{\theta }_0)\) is the limit of \(n_k^{-1}\varvec{\widetilde{H}}_{k}(\varvec{\theta }_0)\) and \(\varvec{\Sigma }^*(\varvec{\beta }_{0})\) is the limit of \(n^{-1}\varvec{H}(\varvec{\beta }_0)\). Since \(\varvec{\Sigma }_k(\varvec{\theta }_0)=\varvec{0}\) for \(k=1,\ldots ,K\) and \(\varvec{\Sigma }_{K+1}(\varvec{\theta }_0)=\varvec{\Sigma }^*(\varvec{\beta }_{0})\), we can conclude from (19) that

The limiting covariance matrix of the distributed estimator \(\hat{\varvec{\beta }}\) is the same as the estimator of \(\varvec{\beta }\) based on the pooled individual-level data Andersen and Gill (1982).

In addition, note that \(\varvec{\Sigma }_{K+1}(\varvec{\theta _0})\), \(\varvec{\Sigma }_{k}(\varvec{\theta _0})\) and \(\varvec{\Sigma }^*(\varvec{\theta _0})\) can be consistently estimated by \(n^{-1}\varvec{\widetilde{H}}( \hat{\varvec{\theta }})\), \(n_k^{-1}\varvec{\widetilde{H}}_{k}(\hat{\varvec{\theta }})\) and \(n_k^{-1}\varvec{I}^*_k(\hat{\varvec{\beta }}_{(k)})\), respectively. In light of (19), we can derive a consistent estimate of the covariance matrix of \(\hat{\varvec{\beta }}\), taking the form

Springer Nature or its licensor (e.g. a society or other partner) holds exclusive rights to this article under a publishing agreement with the author(s) or other rightsholder(s); author self-archiving of the accepted manuscript version of this article is solely governed by the terms of such publishing agreement and applicable law.

Reprints and permissions