1\section{Environment, protocol and assumptions}\label{sec:env}23\subsection{Likelihood, partitions and nested splits}4Let $X_1, \dots, X_N$ be i.i.d.\ with density $f(\cdot;\theta_0)$ with respect to a dominating measure, $\theta_0$ interior to an open $\Theta \subset \R^p$ and $p$ fixed. For a sub-sample $S \subset \{1,\dots,N\}$ write $L_S(\theta) = \abs{S}^{-1}\sum_{i\in S}\log f(X_i;\theta)$, $g_S = \nabla L_S$, $H_S = \nabla^2 L_S$ and let $\thetahat_S$ denote a measurable maximizer of $L_S$ (the MLE on $S$). The centralized MLE is $\thetaN = \thetahat_{\{1,\dots,N\}}$. The index set is partitioned, independently of the data, into $K$ blocks $D_1,\dots,D_K$ of equal size $m = N/K$; \cref{sec:robust} treats unequal sizes. Let $1 = s_1 < s_2 < \dots < s_r$ be integers, each dividing $m$. Level $j$ of worker $k$ consists of the $s_j$ disjoint contiguous sub-blocks of $D_k$ of size $m/s_j$; $T_{s_j,k}$ denotes the average of the $s_j$ MLEs computed on them, so that $T_{1,k} = \thetahat_{D_k}$. The default pattern is $s = (1,2,\dots,2^{r-1})$; the pattern $(1,2,3)$ of \citet{dhaene2015split} is also used. The collection of all sub-blocks at all levels, all blocks and the full sample is the \emph{sub-sample family} $\mathcal S_N$. Asymptotics are along a triangular array $N \to \infty$ with $K = K_N$ and $m \to \infty$.56\subsection{The protocol}\label{sec:protocol}7\Cref{alg:protocol} and \cref{fig:protocol} describe the estimator. Phase~1 is local: worker $k$ fits the $\sum_j s_j$ MLEs of its levels and uploads the $p$-vector8\begin{equation}\label{eq:Jrk}9\Jr{r,k} = \sum_{j=1}^r w_j\, T_{s_j,k},10\end{equation}11where $w = (w_1,\dots,w_r)$ solves the Vandermonde system of \cref{sec:construction}. The coordinator forms the pilot $\pilot{r} = K^{-1}\sum_k \Jr{r,k}$. Each of the $q$ global corrections is one synchronized derivative phase: the coordinator broadcasts the current iterate $\theta$, each worker uploads $g_{D_k}(\theta)$ and the $p(p+1)/2$ distinct entries of $H_{D_k}(\theta)$, the coordinator averages and solves12\begin{equation}\label{eq:newton}13\theta \leftarrow \theta - H_N(\theta)^{-1} g_N(\theta),14\end{equation}15by a linear solve, never by an explicit inverse. We write $\iter{r}{q}$ for the result; $(r,q) = (1,0)$ is the plain average, $(1,1)$ the one-step estimator of \citet{huang2019distributed}, $(2,1)$ the half-block jackknife with one step. Because the derivatives are evaluated at a common data-dependent point, the $q$ derivative phases cannot be merged with phase~1: the protocol has $1 + q$ synchronized exchange phases, $K(1+2q)$ messages, $K(p + q(p + p(p+1)/2))$ uplinked numbers and $Kqp$ downlinked numbers, before headers and serialization (\cref{tab:protocols}). Variants in which the coordinator uses a stale, single-machine or quantized Hessian are analysed in \cref{sec:inexact}.1617\begin{algorithm}[t]18\caption{Order-$r$ local correction with $q$ global Newton corrections}\label{alg:protocol}19\begin{algorithmic}[1]20\Require blocks $D_1,\dots,D_K$ on $K$ workers; splits $s_1 < \dots < s_r$ with weights $w$; rounds $q$21\State \textbf{Phase 1 (workers, in parallel).} For each level $j$, fit the MLE on each of the $s_j$ sub-blocks; set $T_{s_j,k}$ to their average; upload $\Jr{r,k} = \sum_j w_j T_{s_j,k}$ (and a validity flag) \Comment{$p$ numbers up}22\State \textbf{Coordinator.} $\theta \gets K^{-1}\sum_k \Jr{r,k}$23\For{$t = 1,\dots,q$} \Comment{synchronized derivative phase}24 \State broadcast $\theta$ \Comment{$p$ numbers down to each worker}25 \State each worker uploads $g_{D_k}(\theta)$ and the upper triangle of $H_{D_k}(\theta)$ \Comment{$p + p(p+1)/2$ numbers up}26 \State $g \gets K^{-1}\sum_k g_{D_k}(\theta)$, $H \gets K^{-1}\sum_k H_{D_k}(\theta)$; solve $H\,\Delta = g$; $\theta \gets \theta - \Delta$27\EndFor28\State \Return $\iter{r}{q} = \theta$, with the condition number of $-H$ and $\norm{g}$ as diagnostics29\end{algorithmic}30\end{algorithm}3132\begin{figure}[t]33\centering34\begin{tikzpicture}[font=\small, node distance=4mm,35 box/.style={draw, rounded corners=2pt, minimum height=7mm, inner sep=4pt, align=center},36 w/.style={box, fill=blue!6}, c/.style={box, fill=orange!10}, arr/.style={-{Latex[length=2mm]}, thick}]37\node[w] (w1) {worker 1\\ fits on $m, m/s_2,\dots$};38\node[w, right=10mm of w1] (w2) {worker 2};39\node[right=6mm of w2] (dots) {$\cdots$};40\node[w, right=6mm of dots] (wK) {worker $K$};41\node[c, below=9mm of w2, xshift=12mm] (c0) {coordinator: $\theta^{(0)} = K^{-1}\sum_k J_{r,k}$};42\node[c, below=7mm of c0] (c1) {broadcast $\theta$; aggregate $g_{D_k}(\theta), H_{D_k}(\theta)$; solve; repeat $q$ times};43\node[c, below=7mm of c1] (c2) {$\theta^{(q)}_r$};44\foreach \x in {w1,w2,wK} {\draw[arr] (\x.south) -- node[pos=.5, left=1pt, font=\scriptsize]{} (c0.north -| \x.south);}45\draw[arr] (c0) -- node[right=2pt, font=\scriptsize]{$q$ synchronized derivative phases: $p$ down, $p + p(p+1)/2$ up, per worker and phase} (c1);46\draw[arr] (c1) -- (c2);47\node[below=1mm of c2, font=\scriptsize, text width=12.5cm, align=center] {Phase 1 uplink: $p$ numbers per worker (plus $p(p+1)/2$ if the local Hessian is shipped for the stale-curvature variant). Arrows are logical exchanges, not measured latency; \cref{sec:bench} reports measured times and bytes.};48\end{tikzpicture}49\caption{Coordinator--worker protocol with exact uplink and downlink accounting.}\label{fig:protocol}50\end{figure}5152\subsection{Assumptions}\label{sec:assumptions}53The theorems are conditional on the following assumptions, stated for the sub-sample family $\mathcal S_N$. They are high-level: they are implied by standard primitive conditions for regular parametric models (smoothness, moment and identification conditions of the type used by \citealp{zhang2013communication} and \citealp{huang2019distributed}, together with the validity of the formal bias expansion, \citealp{bhattacharya1978validity}), and they are verified exactly in the exponential model of \cref{sec:exact}.5455\begin{assumption}[Uniform stochastic expansion]\label{ass:lin}56There exist a measurable $\psi$ with $\E\psi(X) = 0$ and $\Sigma_\psi := \Var\psi(X)$ finite and positive definite, an integer $n_0$ and a constant $C < \infty$ such that for every $S \in \mathcal S_N$ with $\abs{S} \ge n_0$,57\begin{equation}\label{eq:lin}58\thetahat_S - \theta_0 = \abs{S}^{-1}\sum_{i\in S}\psi(X_i) + R_S, \qquad \E\norm{R_S}^2 \le C\,\abs{S}^{-2}.59\end{equation}60\end{assumption}61For a correctly specified regular likelihood $\psi = I(\theta_0)^{-1}\nabla_\theta\log f(X;\theta_0)$ and $\Sigma_\psi = I(\theta_0)^{-1}$; for a misspecified model $\theta_0$ is the pseudo-true value, $\psi = -\bar H(\theta_0)^{-1}\nabla_\theta\log f(X;\theta_0)$ and $\Sigma_\psi$ is the sandwich (\cref{sec:robust}). The same $\psi$ serves every nested sub-sample; this is the only property of the influence function used.6263\begin{assumption}[Bias expansion of order $r$]\label{ass:bias}64There exist vectors $b_1,\dots,b_{r-1}$ and $C_b < \infty$ with65\begin{equation}\label{eq:bias}66\E R_S = \sum_{a=1}^{r-1} b_a\,\abs{S}^{-a} + \rho_{\abs{S}}, \qquad \norm{\rho_n} \le C_b\,n^{-r}, \qquad S \in \mathcal S_N,\ \abs{S}\ge n_0.67\end{equation}68\end{assumption}69For $r = 1$ the assumption only asks $\norm{\E R_S} \le C_b\abs{S}^{-1}$, which follows from \cref{ass:lin} by Jensen's inequality. The existence of the expansion for $r \ge 2$ is not implied by consistency or asymptotic normality; it requires smoothness of the likelihood to order $r+1$ and moment conditions, and it fails for estimators with non-smooth criteria.7071\begin{assumption}[Local Newton regularity]\label{ass:newton}72With probability tending to one, $\thetaN$ is an interior zero of $g_N$; there are a fixed convex neighbourhood $U$ of $\theta_0$ and random variables $M_N, \Lambda_N = \op(1)$ such that $\sup_{\theta\in U}\norm{H_N(\theta)^{-1}} \le M_N$ and $\norm{H_N(u) - H_N(v)} \le \Lambda_N\norm{u - v}$ for $u, v \in U$. Moreover $\sqrt{N}(\thetaN - \theta_0) \dto \mathcal N(0,\Sigma_\psi)$.73\end{assumption}7475\begin{assumption}[Observation-level curvature smoothness; used only in \cref{thm:inexact}(a)]\label{ass:third}76$h(x;\theta) := \nabla^2_\theta\log f(x;\theta)$ is continuously differentiable on $U$ with $\norm{\nabla h(x;\theta) - \nabla h(x;\theta')} \le L(x)\norm{\theta - \theta'}$, $\E L(X)^2 < \infty$, $\E\norm{h(X;\theta_0)}^2 < \infty$, $\E\norm{\nabla h(X;\theta_0)}^2 < \infty$, and $\E\norm{\thetahat_{D_k} - \theta_0}^4 \le C_4 m^{-2}$.77\end{assumption}7879\begin{remark}80The assumptions fix $p$, require independent identically distributed blocks and well-behaved local optimizers. Constants such as $C$, $M_N$ and $\Lambda_N$ depend on $p$ and on the parametrization; none of the exponents below should be read as valid when $p$ grows with $N$. The condition $m \to \infty$ is essential: all statements are about a triangular array in which the number of partitions and the local sample size grow together, not about a fixed data set with a maximal number of machines.81\end{remark}82