A Kernel Perspective on Distillation-based Collaborative Learning

Authors: sejun  park
kihun  hong
ganguk  hwang
Venue: NeurIPS_cc_2024_Conference
Published: 2024-05-14 08:24:51

Abstract: Over the past decade, there is a growing interest in collaborative learning that can enhance AI models of multiple parties.
However, it is still challenging to enhance performance them without sharing private data and models from individual parties.
One recent promising approach is to develop distillation-based algorithms that exploit unlabeled public data but the results are still unsatisfactory in both theory and practice.
To tackle this problem, we rigorously analyze a representative distillation-based algorithm in the view of kernel regression.
This work provides the first theoretical results to prove the (nearly) minimax optimality of the nonparametric collaborative learning algorithm that does not directly share local data or models in massively distributed statistically heterogeneous environments.
Inspired by our theoretical results, we also propose a practical distillation-based collaborative learning algorithm based on neural network architecture.
Our algorithm successfully bridges the gap between our theoretical assumptions and practical settings with neural networks through feature kernel matching.
We simulate various regression tasks to verify our theory and demonstrate the practical feasibility of our proposed algorithm.

Main paper: 
 introduction: 
 Collaborative learning of AI models in decentralized settings is an important problem covered in various fields of machine learning such as distributed learning [10, 70] , Federated Learning (FL) [23] , peer-to-peer learning [3] , and miscellaneous collaborative learning [47] . In particular, this theme has been most actively discussed in the context of FL [24, 31, 45, 61] . In this context, each local party is typically viewed as a subordinate entity within the collective learning system. For example, most FL algorithms mandate the exchange of local AI model information among participating local parties. Under this scheme, local AI models are usually subjected to restrictions in their architecture. However, from the perspective of collaboration, each local party may have to be regarded as an independent learning agent, meaning they are not obligated to fully share their model information.In short, the model (or parameter) exchange in FL algorithms can emerge as a critical issue in collaborative learning.Fundamentally, addressing this issue necessitates an alternative medium for sharing learning information distinct from model exchange. Indeed, Distillation-based Collaborative Learning (DCL) [14, 30, 44] provides a good answer. In these algorithms, local training information is shared via the outcomes of AI models on additional unlabeled public data. The collected information is then utilized for knowledge distillation [21] to each local AI model. As mentioned in [14, 48] , this procedure is agnostic to model heterogeneity and avoids the direct sharing of local AI model information. This is a key advantage that distinguishes DCL from traditional FL.Despite its pioneering nature and potential utility, DCL has not been sufficiently explored. A significant reason for this is the lack of theoretical understanding regarding knowledge distillation and its effectiveness in massively distributed statistically heterogeneous environments. Our work stems from the fundamental question of whether DCL algorithms can be theoretically effective in these settings. Inspired by [48, 58] , we analyze FedMD [30, 41] , the most standard DCL algorithm from a nonparametric perspective. Specifically, we adopt an operator-theoretic approach [5, 15, 37, 53, 65] to obtain an upper rate of convergence for the nonparametric version of FedMD (named DCL-KR) in the expected sense. Remarkably, our analysis reveals that DCL-KR achieves a nearly minimax optimal convergence rate, where the prefactor is independent of the number of participating local parties. It is worth noting that DCL-KR is the first nearly minimax optimal collaborative learning algorithm that does not directly share local data or models in massively distributed statistically heterogeneous environments. The novelty of our theoretical results and their comparison to prior works are provided in Section 2 and 3.Nevertheless, our theoretical analysis does not fully demonstrate the efficacy of DCL algorithms based on neural network architectures. Instead, our theoretical results serve as inspiration for designing a novel DCL algorithm for regression that refines existing approaches. Consequently, we propose a Distillation-based Collaborative Learning algorithm over heterogeneous Neural Networks (named DCL-NN) for regression tasks. DCL-NN leverages kernel matching to align the feature kernels from the last hidden layer of each local AI model with an ensemble kernel. This procedure brings heterogeneous neural networks into the regime of DCL-KR.Finally, we conduct experiments on DCL-KR and DCL-NN. To illustrate the superiority of our algorithms, we compare them with several baselines on various regression tasks. Experimental results show that DCL-KR achieves the same performance as the centralized model, even beyond the theoretical results. We also observe that DCL-NN significantly outperforms previous DCL frameworks in most settings.In summary, our contributions are as follows:1. In Section 3, we theoretically prove that a nonparametric version of the most standard distillationbased collaborative learning algorithm (named DCL-KR) is nearly minimax optimal in massively distributed statistically heterogeneous environments. 2. Inspired by the results provided in Section 3, we propose a distillation-based collaborative learning algorithm with heterogeneous neural networks (named DCL-NN) in Section 4. 3. In Section 5, we conduct experiments to empirically confirm our theoretical results and show the practical feasibility of our proposed algorithms. 
related work: 
 Federated Learning Most FL algorithms [45] communicate model parameters for collaboration. This approach has been extensively studied under various constraints, including data privacy [1] , statistical heterogeneity [24, 31] , communication efficiency [51] , personalization [13, 59] , and robustness [25] . While it has been successful both theoretically and experimentally, this type of FL is limited in terms of the privacy and flexibility of local AI models, as the algorithms directly access the structures and parameters of the local models. Our study focuses on distillation-based collaborative learning, where the privacy and flexibility of local AI models are fully guaranteed.The type of algorithms we investigate operates by communicating the functional information of local AI models. These algorithms typically assume the availability of additional public data points. In this case, the outcomes of local models on the public dataset are used for collaboration. For instance, Li and Wang [30] , Lin et al. [41] , Park et al. [48] iteratively collect predictions of local models on the public dataset and then aggregate them into a naive ensemble (with or without a fixed linear transformation) to distribute. On the other hand, Cho et al. [7] , Zhang et al. [69] , Fan et al. [14] apply personalized ensemble strategies by additionally learning the mutual trust between models. Makhija et al. [44] propose FedHeNN, which distills training information in the form of matching feature kernels instead of the predictions of local AI models on the public data. Both FedHeNN and DCL-NN utilize centered kernel alignment [8] to match feature kernels of local models, but DCL-NN uses the ensemble distillation for predictions as well. Thus, DCL-NN enables parties to learn from the entire input space.Table 1 : Comparative analysis of decentralized environments for (nearly) minimax optimality of representative collaborative learning algorithms with kernel regression. nFedAvg indicates the nonparametric version of FedAvg in [58] . Note that IED [48] achieves a weaker version of minimax optimality.interaction local data massively non-i.i.d.& Methods method privacy distributed unbalanced DKRR [38] divide-and-conquer DC-NY [66] divide-and-conquer DKRR-CM [40] model exchange DKRR-RF-CM [43] model exchange DKRR-NY-CM [67] model exchange nFedAvg [58] model exchange IED* [48] knowledge distillation DCL-KR (ours) knowledge distillationDecentralized Learning with Kernel Regression A number of studies have investigated the minimax optimal rate of regularized kernel regression algorithms such as kernel ridge regression and gradient descent-based kernel regression with early stopping [5, 15, 37, 65] . In particular, over the past decade, the growing interest in decentralized learning has led to active research in the generalization analysis of decentralized kernel regression. While divide-and-conquer algorithms [34, 38, 66, 70] play a significant role in this research flow, most of them fail to account for statistical heterogeneity and massively distributed cases, along with privacy preservation, which has received a lot of attention recently. On the other hand, decentralized kernel regression algorithms with multiple communication rounds [40, 43, 48, 58, 67] achieve superior theoretical results compared to the divide-and-conquer algorithms. However, the discussions of these algorithms primarily focus on the efficiency of resource costs [40, 43, 67] , while research on relaxing environmental constraints has been scarce. For example, most of these works assume a limited number of parties to prove the optimality in a minimax sense.To the best of our knowledge, [48, 58] stand as the only investigations that consider general decentralized environments. Similar to our work, Park et al. [48] study the convergence rate of distillation-based collaborative learning with kernel regression. However, their results demonstrate a weaker version of minimax optimality and do not cover statistically heterogeneous environments. In this regard, Su et al. [58] offer a promising methodology. They analyze nonparametric versions of FedAvg [45] and FedProx [31] , representative FL algorithms involving model exchange, in general decentralized environments such as statistically heterogeneous and massively distributed scenarios. In this work, we extend their methodology to analyze FedMD [30, 41] from a nonparametric perspective in massively distributed statistically heterogeneous environments. We summarize the comparison between our work and prior studies in Table 1 . Note that algorithms that do not employ Nyström scheme (including nonparametric FedAvg [58] ) fail to preserve local data privacy due to the inherent characteristics of kernel regression. On the other hand, DC-NY [66] and DKRR-NY-CM [67] can achieve the local data privacy preservation by utilizing the public data as Nyström centers. 
dcl-kr: a nonparametric view of fedmd: 
 In this section, we establish the theory of a nonparametric version of FedMD [30, 41] , the most standard distillation-based collaborative learning algorithm.Let ρ x,y = ρ x • ρ y|x be a Borel probability measure on X × R where X is a compact subset of R d and we assume the support of ρ x is X . The goal of the regression problem is to find a minimizer of the population risk, i.e.,min h:X →R E(h), E(h) := 1 2 E (x,y)∼ρx,y |y -h(x)| 2 .Then, the function f * 0 : X → R defined byx 0 → E y∼ρ y|x (•|x0) [y], x 0 ∈ X is a target function.Let k : X × X → R be a Mercer kernel [9] where κ := (sup x∈X k(x, x)) 1/2 < ∞ and H k be a reproducing kernel Hilbert space associated to k. We set k x := k(•, x) and the covariance operator T k,ν : H k → H k with respect to any Borel probability measure ν on X defined asT k,ν h = X h(x)k x dν(x).Then we can see that T k,ν = ι ν ι ν where ι ν :H k → L 2 ν is a natural embedding, L 2 ν = L 2 (X , ν) denotes the L 2 space, and a superscript denotes the adjoint operator of a given operator. We also define the sampling operator S D :H k → R n by h → [h(x 1 ), • • • , h(x n )] and T k,X := S D S D when D = {(x 1 , y 1 ), • • • , (x n , y n )} with X = {x 1 , • • • , x n }is given. Since S D depends only on data inputs X, we can define the sampling operator for unlabeled datasets in the same way. See Appendix A.1 for further details.Given a dataset D = {(x 1 , y 1 ), • • • , (x n , y n )} generated from ρ x,y , consider the empirical risk E D : H k → R given by E D (h) = 1 2 S D h -y 2 2 where y = [y 1 , • • • , y n ] . Here, • 2 denotes a scaled Euclidean norm v 2 = ( 1 n n i=1 v 2 i ) 1/2 . From the functional derivative ∇ E D (h) = S D (S D h -y), the gradient descent scheme becomes ν 1 = 0, ν t+1 = ν t -η t S D (S D ν t -y) (t = 1, 2, • • • )where {η t } t∈N is a set of learning rates. In this work, we set η t = η, t ∈ N for a fixed η ∈ (0, 1/κ 2 ). Then, a simple calculation gives ν t → S D (S D S D ) -1 y as t → ∞ provided that the operator S D S D is invertible. The limit is known as the minimum norm interpolation [49] of D. Since the interpolation regressor generalizes poorly unless there is no noise [32, 39] , early stopping strategies are usually applied to avoid the overfitting issue. With adequate stopping rules, gradient descent-based kernel regression has an optimal rate in a minimax sense [36, 37, 65] .From now on, we consider the setting that there are m parties and the ith party has a private local dataD i = {(x j i , y j i ) : j = 1, • • • , n i } for i = 1, • • • , m.Assume that all data D = m i=1 D i are i.i.d. with the distribution ρ x,y but each local dataset does not need to have the same distribution. Let Z = {z 1 , • • • , z n0 } ⊂ X be the additional public inputs. The goal of all parties is to have their models that perform well on the distribution ρ x . In other words, each party expects to be able to make good predictions not only for its local data distribution but also for unseen data distribution through collaborative learning.Similar to [58] , we construct a nonparametric version of FedMD (called DCL-KR), which is presented in Algorithm 1. In Algorithm 1, G i is a one-step local gradient descent update on E Di , i.e.,G i h = h -ηS Di (S Di h -y i ) where y i = [y 1 i , • • • , y ni i ] . Similarly, Gt is a one-step gradient descent update on E (Z,yp,t) , i.e., Gt h = h -ηS Z (S Z h -y p,t ).In this subsection, we show the nearly minimax optimality of DCL-KR. To derive theoretical results, we assume the following conditions regarding regularity of noise, the kernel k, and the target function f * 0 as below. Assumption 3.1. We assume E y∼ρy y 2 < ∞ andexp |y -f * 0 (x)| M - |y -f * 0 (x)| M -1 dρ y|x (y|x) ≤ γ 2 2M 2 , ∀x ∈ Xwhere M and γ are positive constants. i,0 = 0 for i = 1, • • • , m. 3: for t = 0, • • • , T -1 do 4: for party i = 1, • • • , m do 5:Update the local model E times by gradient descent on the empirical risk E Dif i,t ← G E i f i,t . 6:Upload the local predictions on Z to the servery i p,t = S Z f i,t . 7:end for 8:The server aggregates the local predictions to compute the consensus predictiony p,t = m i=1 n i n y i p,tand then distributes y p,t to all local parties. 9:For party i(i = 1, • • • , m), update the local model by infinitely many iterations of gradient descent on the empirical risk E (Z,yp,t)EQUATIONwith an initialization g i,t chosen from a subspace spanned byk z 1 , • • • , k z n 0 . 10: end for Assumption 3.2. Let λ 1 ≥ λ 2 ≥ • • • > 0 be eigenvalues of T k,ρx .There are fixed positive constants C s and c s such thatc s i -1/s ≤ λ i ≤ C s i -1/s , ∀i ∈ Nfor some s ∈ (0, 1).Assumption 3.3. The target function f * 0 satisfiesf * 0 ∈ h ∈ H k : h = T r-1/2 k,ρx g where g H k ≤ R for some r ∈ [ 1 2 , 1] where T r-1/2 k,ρxis the (r -1/2) power of operator T k,ρx and R > 0 is a fixed constant. In particular, f * 0 ∈ H k .The above assumptions determine the minimax lower rate [5] and are standard assumptions in many prior works [5, 15, 33, 35] . In detail,• Assumption 3.1 implies that the noise is not excessively large. This assumption is a general noise condition that encompasses a wide range of cases. For instance, noise with Bernstein condition such as sub-Gaussian noise satisfies Assumption 3.1.• Assumption 3.2 is about the eigenvalue decay of T k,ρx . From this assumption, one can derive bounds on the effective dimension that is related to covering and entropy number conditions [15] .• Assumption 3.3 is related to the regularity of the target function, specifically how well the RKHS induced by the kernel k represents the target function.Under these assumptions, we can theoretically show the performance guarantee of DCL-KR. The proof is provided in Appendix A.2. Note thatE(h) -E(f * 0 ) = 1 2 ι ρx (h -f * 0 ) 2 L 2 ρxis the excess risk of a regressor h and so the quantity ι ρx (h -f * 0 ) L 2 ρx indicates the generalization ability of h. EQUATIONDCL-KR gives the performance guaranteeE ι ρx (f i,T -f * 0 ) L 2 ρx ≤ C • B r n -r 2r+s log n for all i = 1, • • • , m where η ∈ (0, 1/κ 2) is a fixed learning rate, T is an adequate stopping rule, and the prefactor C does not depend on B, m, and n.Since the convergence rate n -r 2r+s is the minimax lower rate under Assumption 3.1, 3.2, and 3.3, Theorem 3.4 implies that DCL-KR has an almost same convergence rate as the minimax optimal central training when there are sufficiently many public inputs. To the best of our knowledge, this is the first work to prove the (nearly) minimax optimality of a collaborative learning algorithm that does not directly share local data or models in massively distributed statistically heterogeneous environments. For example, divide-and-conquer algorithms work for limited m. Specifically, DC-NY [66] assumes m ≤ O(n 2r-1 2r+s ) and DKRR-NY-CM [67] assumes m ≤ O(n 2r+s-1 2r+s ). However, Theorem 3.4 does not require any condition on m. Moreover, Theorem 3.4 deals with a more general setting than the theory in [48, 58] . For example, Su et al. [58] only cover r = 1 2 of Assumption 3.3. On the other hand, Park et al. [48] do not consider Assumption 3.2 which gives a finer result. Compared with [48] , we also reduce the required size of public inputs and drop the statistical homogeneity condition.The convergence rate in Theorem 3.4 has an additional factor log n compared with a minimax lower rate [5, 15] , but this logarithm term grows slower than any polynomial. Note that an additional logarithm term commonly appears in the context of gradient descent-based kernel regression with Nyström scheme [35, 36] . Theorem 3.4 allows that the public input distribution ρx can be different from the local input distribution ρ x . It is natural that the condition (2) is required since ρx should cover ρ x for fully distilling training information. We can see that the discrepancy between ρ x and ρx affects the upper bound in Theorem 3.4 as the multiplication of B r . We can remove B r in the upper bound by increasing public inputs. See Appendix A.3 for details.In the proof of Theorem 3.4, we decompose the term ι ρx (f i,T -f * 0 ) into four parts, say (I), (II), (III), and (IV) (see Eq. ( 7)). The proof is to bound the norms of these terms. Note that DCL-KR can also be understood as a Nyström version of nonparametric FedAvg [58] from the recurrence relation (6) .(I) and (II) appear similarly in [58] , except that (I) and (II) incorporate projections. To handle these terms, we reinterpret the proof presented in [58] in operator form instead of matrix form and extend it to our setting. We obtain a norm bound of (II) containing a quantity linked to the local Rademacher complexity. (Appendix A.2.2 and A.2.3)Comparing with [58] , (III) and (IV) are additional terms induced by the procedure that distills functional information from the local regressors. We apply techniques used in [35, 48, 53] to bound (III) and (IV). (Appendix A.2.4) Note that previous works applying local Rademacher complexity-based stopping rule [50, 58] deal with the case of r = 1 2 only. In this work, we set a new stopping rule T which is an extension of previous works [50, 58] and prove an extended version (Lemma A.6) of a well-known property [60] . As a result, our theory covers r ∈ [ 1 2 , 1] which affects the minimax lower rate. (Appendix A.2.5) 
dcl-nn algorithm: 
 In this section, we retain the problem setting from Section 3 but employ heterogeneous neural networks as the local models. Based on the theoretical results in Section 3, we propose a novel distillation-based collaborative learning algorithm DCL-NN across heterogeneous neural networks in a decentralized setting.A key factor contributing to the successful theoretical guarantee of DCL-KR lies not only in the linearity of kernel regression but also in the equality of kernels across local models. In fact, the public data predictions can vary in different directions, even if the same training data points are used when kernels differ (See Appendix B). Therefore, we match the kernels of local AI models. Specifically, we use linear feature kernels [18, 64] induced by the features from the last hidden layers of local AI models for kernel matching. For example, for a neural network f :X → R where f (•) = w g(•) + b, g : X → R c , w ∈ R c, and b ∈ R we useEQUATIONas the feature kernel of f . Through this idea, we can bring the setting closer to the regime of DCL-KR. Note that our theoretical results suggest that the target kernel should be a good kernel. Indeed, we observe that the naive ensemble For this purpose, we introduce Centered Kernel Alignment (CKA) [8] as a kernel similarity measure. CKA is a typical measure associated with the similarity of two representations of neural networks [27] and is often used for kernel matching in neural networks [44] . To compute empirical CKA between two kernels k 1 and k 2 on inputs {c 1 , • • • , c p }, we first calculate the Gram matricesEQUATIONK 1 = [k 1 (c j1 , c j2 )] 1≤j1,j2≤p and K 2 = [k 2 (c j1 , c j2 )] 1≤j1,j2≤p. We then compute the empirical CKA viaCKA(k 1 , k 2 ) = HSIC(K 1 , K 2 ) HSIC(K 1 , K 1 ) HSIC(K 2 , K 2 ).Here, HSIC is an estimator of the Hilbert-Schmidt Independence Criterion (HSIC) defined asHSIC(K 1 , K 2 ) = 1 (p -1) 2 tr(K 1 HK 2 H)where H := I p -1 p 11 is the centering matrix. In the kernel distillation procedure, the ith local party maximizes CKA(k fi , k) on public inputs Z (i = 1, • • • , m). Notably, this procedure requires only a single communication round for exchanging pairwise feature kernel values on public inputs, ensuring that our algorithm operates exclusively within the function space.After the kernel distillation procedure, all local AI models have similar feature kernels up to constants. So we can follow an analogous process as in DCL-KR. Note that we perform learning rate scaling described in Appendix B to compensate the kernel scale difference. It makes the impact of local iterations consistent. We also provide the complete algorithm (Algorithm 2) and further details for Section 4 in Appendix B. 
experiments: 
 In this section, we evaluate the performance of DCL-KR and DCL-NN. We compare them with baselines on various regression tasks.Datasets We use the following six regression datasets to evaluate the performance. Target variables are one-dimensional in all datasets. (1) Toy-1D [33] and (2) Toy-3D [6] are synthetic datasets with one-dimensional and three-dimensional inputs, respectively. (3) Energy is a tabular dataset from the UCI database [12] to predict appliances energy use with 28 features. (4) RotatedMNIST is an image dataset where it aims to predict the rotation angles for given rotated images of the MNIST [11] images. (5) UTKFace [71] and (6) IMDB-WIKI [42, 52] are image datasets for age estimation. We compare kernel machine-based collaborative learning algorithms on two datasets Toy-1D and Toy-3D. On the other hand, we compare neural network-based collaborative learning algorithms on five datasets Toy-3D, Energy, RotatedMNIST, UTKFace, and IMDB-WIKI.Baselines We compare DCL-KR with two central kernel regression models to verify our theoretical results. These two central models have the minimax optimal convergence rate. We also utilize existing decentralized kernel regression algorithms that does not directly share local data and models (DC-NY [66] , DKRR-NY-CM [67] , IED [48] ) as baselines for DCL-KR. On the other hand, we adopt FedMD with unlabeled public inputs [30, 41] , FedHeNN [44] , and KT-pFL [69] as baselines for DCL-NN.Setup The number of parties ranges from 10 to 100 for kernel machine-based algorithms and is 50 for neural network-based algorithms. We construct statistically heterogeneous decentralized environments with Algorithm 3. For neural network-based algorithms, we use 4 different neural network architectures for local models in all settings. For instance, we use ResNet-18, ResNet-34, ResNet-50 [20] , and MobileNetv2 [54] for large-scale image datasets. We utilize the average of Root Mean Squared Errors (RMSEs) of the local AI models on a test dataset as a performance metric. The test data points have the same distribution as the whole local data distribution. We apply FedMD with a few communication rounds for pretraining of DCL-NN. See Appendix C for detailed experimental configurations.The performance of DCL-KR and its baselines is presented in Figure 1 . We set the number of parties m = 10, 20, • • • , 100, the number of private data points n = 50m, and the number of public inputs n 0 = n 1 2r+s (log 10 n) 3 . We first set ρ x = ρx , i.e., the public data distribution is the same as the entire local input distribution. As shown in Figure 1 , DCL-KR outperforms the baselines in all experimental settings and achieves comparable performance to the central models. This result implies that DCL-KR has not only the nearly optimal convergence rate but also the same performance as central kernel regression models. In contrast, DC-NY and DKRR-NY-CM exhibit significantly lower performance compared with DCL-KR in massively distributed environments where their theory does not cover. IED does not show a significant performance drop in massively distributed environments even though its theory is built on the statistical homogeneity condition of local data distributions.To further compare the performance of DCL-KR and IED, which show similar results to central models, we analyze the effect of n 0 and ρx on their performance (Figure 2 and 3 ). Figure 2 illustrates that, as expected from the theoretical results, IED requires more public inputs than DCL-KR to achieve good performance. Moreover, when there is a public distribution shift, DCL-KR maintains its convergence rate, whereas the convergence rate of IED deteriorates. (See Appendix C.3.3 for experimental details.) Overall, our experiments validate the theoretical results of DCL-KR and demonstsrate its superiority over previous results. For additional experimental results and analyses, please refer to Appendix C.3.Table 2 shows the performance of DCL-NN and baselines on five regression tasks. We also present the performance of standalone models and centralized models to assess the performance of the collaborative algorithms. For some cases exhibiting training instability, we report the best test error (marked with asterisks) observed across all communication rounds, while relying on a fixed number of communication rounds for the other cases.As can be seen in Table 2 , DCL-NN outperforms the baselines on all regression tasks. Note that FedHeNN employs kernel matching similar to DCL-NN, but it lacks supervision of label prediction through collaboration, resulting in insufficient performance improvement compared to standalone models. Given the superior performance of DCL-NN, it is evident that incorporating supervised learning for label prediction alongside kernel matching is desirable. On the other hand, while FedMD performs significantly better than standalone models, the performance of DCL-NN is consistently better. Considering that we utilize FedMD for pretraining of DCL-NN, we can see that it performs better than FedMD-only collaborative learning by first training local models with FedMD and then using DCL-NN. In conclusion, the experimental results support the practical effectiveness and superiority of DCL-NN over baselines.Kernel Distillation Procedure To verify the necessity of kernel distillation, we examine the changes in the performance of local feature kernels and the CKA between them during the kernel distillation procedure. We conduct this experiment on UTKFace. We utilize the RMSE of a kernel linear regression model trained on all local data as a kernel performance measure. The results are presented in Figure 4 . As shown in Figure 4 , both kernel performance and CKA undergo a temporary degradation due to the change of the objective function at the initial stages. However, as training progresses, both metrics recover and kernel performance surpasses its initial level. Since kernel distillation aims to ensure that all local feature kernels are similar with high performance, the experimental results verify the effectiveness of kernel distillation.Table 2 : Performance comparison of FedMD, FedHeNN, KT-pFL, and DCL-NN on five datasets. The values are presented as the average of RMSEs along with standard deviations. For calibration, the performance of standalone models and centralized models is also provided.  
conclusions: 
 In this work, we analyze distillation-based collaborative learning from a nonparametric perspective and propose DCL-NN, a practical algorithm as an extension. We demonstrate that DCL-KR, a nonparametric version of FedMD, has a nearly minimax optimal convergence rate in massively distributed statistically heterogeneous environments. Inspired by DCL-KR, we propose DCL-NN, a novel distillation-based collaborative learning algorithm for heterogeneous neural networks. Our experiments confirm the theoretical results of DCL-KR and demonstrate the practical effectiveness of DCL-NN. For a discussion of the limitations of our work, please refer to Appendix D.Broader Impact Our work explores the methodologies of collaborative learning under data and model privacy preservation. In this regard, our research holds the potential to positively impact the facilitation of collaboration among AI models without raising concerns about information disclosure. On the other hand, our work does not pose any particularly noteworthy negative consequences, given its aim to contribute to the advancement of the general field of machine learning. d-dimensional Euclidean space X the input space contained in R d C(X )the collection of all continuous functions from X into R ρ x,y the data generating distribution on X × R ρ x , ρ y the marginal distribution of ρ x,y on X and R, respectivelyρ y|x (•|x 0 ) the conditional distribution on R w.r.t. x 0 ∈ X and ρ x,y ρx the public input distribution on X k a given Mercer kernel k x k(•, x) κ (sup x∈X k(x, x)) 1/2 L 2 ν the L 2 space on X w.r.t. measure ν H k a reproducing kernel Hilbert space associated to kernel k T k,ν the covariance operator on H k w.r.t. measure ν, T k,ν : h → X h(x)k x dν(x) ι ν a natural embedding from H k into L 2 ν S D a sampling operator from H k into R n , S D : h → [h(x 1 ), • • • , h(x n )] where D = {(x i , y i )} n i=1 S Z a sampling operator from H k into R n , S Z : h → [h(z 1 ), • • • , h(z n )] where Z = {z i } n i=1 T k,X S X S X E(h) the population risk of h, E(h) = 1 2 E (x,y)∼ρx,y |y -h(x)| 2 E D (h) the empirical risk of h over D = {(x i , y i )} n i=1 , E D (h) = 1 2 S D h -y 2 2 f * 0 a target function from X into R defined as f * 0 (x 0 ) = E y∼ρ y|x (•|x0) [y], x 0 ∈ X η a learning rate, η ∈ (0, 1/κ 2 ) m the number of parties D i private local data of the ith party, {(x j i , y j i ) : j = 1, • • • , n i } (i = 1, • • • , m) X i inputs of D i , {x j i : j = 1, • • • , n i } (i = 1, • • • , m) y i labels of D i , [y 1 i , • • • , y ni i ] (i = 1, • • • , m) D D = m i=1 D i Z unlabeled public data, {z 1 , • • • , z n0 } M, γthe parameters related to the regularity of noise (Assumption 3.1) (5) . s, C s , c s the parameters related to the eigenvalue decay of T k,ρx (Assumption 3.2){(λ i , φ i )} ∞ i=1 eigenvalues and eigenvectors of T k,ρx such that λ 1 ≥ λ 2 ≥ • • • > 0 from Mercer's representationC s C s s /(1 -s) N ν (λ) tr(T k,ν (T k,ν + λI) -1 ) where ν is a probability measure N (λ) N ρx (λ) r, Rthe parameters related to the regularity of f * 0 (Assumption 3. A Details on Section 3Before we start the proof of Theorem 3.4, we present basic notions. 
experiments compute resources: 
 Question: For each experiment, does the paper provide sufficient information on the computer resources (type of compute workers, memory, time of execution) needed to reproduce the experiments?Answer: [Yes] Justification: Descriptions of computer resources are provided in Appendix C.2. 
code of ethics: 
 Question: Does the research conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics https://neurips.cc/public/EthicsGuidelines? 
broader impacts: 
 Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?Answer: [Yes] Justification: See Section 6. 
 
 Appendix: In Subsection 3.1, the reproducing kernel Hilbert space H k is a subset of C(X ), i.e., all elements in H k are continuous [57] . Sinceι ρx h(•) = ι ρx h, k • H k = h, ι ρx k • L 2 ρx = X h(x)k(•, x) dρ x (x),we have T k,ρx = ι ρx ι ρx . The compactness of ι ρx [57] gives the fact that T k,ρx is compact, self-adjoint, and positive. Furthermore, Mercer's theorem [62] gives a Mercer representationEQUATIONThe fact that H k ⊂ C(X ) and T k,ρx = ι ρx ι ρx implies the injectivity of T k,ρx and so λ i = 0 for all i ∈ N. We define N ν (λ) := tr(T k,ν (T k,ν + λI) -1 ) for any probability measure ν. For convenience, N (λ) := N ρx (λ). From [5, 15] , we haveN (λ) ≤ C s λ -s where C s := C s s /(1 -s) and N ν (λ) ≤ κ 2 λ -1 . Given a dataset D = {(x i , y i )} n i=1, a similar argument as above givesS D : c = [c 1 , • • • , c n ] → 1 n n i=1 c i k x i and T k,X : h → 1 n n i=1 h(x i )k x i . Note that Assumption 3.1 implies that E(f * 0 ) = 1 2 E (x,y)∼ρx,y |y -f * 0 (x)| 2 ≤ γ 2 2 < ∞.We haveE (x,y)∼ρx,y |y -h(x)| 2 = E x∼ρx |h(x) -f * 0 (x)| 2 + E (x,y)∼ρx,y |y -f * 0 (x)| 2and so the excess risk becomesE(h) -E(f * 0 ) = 1 2 E x∼ρx |h(x) -f * 0 (x)| 2 = 1 2 ι ρx (h -f * 0 ) 2 L 2 ρx . Therefore, ι ρx (h -f * 0 ) 2 L 2 ρxindicates the generalization ability of h.Table 3 presents meaning of some notations.A.2 Proof of Theorem 3.4Without loss of generality, we assume n ∧ n 0 ≥ κ 2 e.Consider a subspace W of H k spanned by {k z 1 , • • • , k z n 0 }. We first show that for a fixed h * ∈ H k and a gradient update Gu = u -ηS Z (S Z u -S Z h * ) we have G t u 1 → P Z h * as t → ∞ for any u 1 ∈ W where P Z is an orthogonal projection onto the subspace W . Set u t+1 = Gu t for t ≥ 1.Then ut+1 = (I -ηS Z SZ )ut + ηS Z SZ h * = (I -ηS Z SZ ) t u1 + t-1 k=0 (I -ηS Z SZ ) k ηS Z SZ h * . Since S Z h * = S Z P Z h * , we have t-1 k=0 (I -ηS Z SZ ) k ηS Z SZ h * = t-1 k=0 (I -ηS Z SZ ) k ηS Z SZ PZ h * = PZ h * -(I -ηS Z SZ ) t PZ h * .Note that there exists{z 1 , • • • , zñ0 } ⊂ Z such that {k z1 , • • • , k zñ 0 } is a basis of W . Define a matrix B =    b11 • • • b1ñ 0 . . . . . . . . . bn 0 1 • • • bn 0 ñ0    ∈ R n 0 ×ñ 0 such that k z i = ñ0 j=1 b ij k zj . Then K Z Z = BK Z Z where K Z Z =    k(z 1 , z1 ) • • • k(z 1 , zñ 0 ) . . . . . . . . . k(z n 0 , z1 ) • • • k(z n 0 , zñ 0 )    ∈ R n 0 ×ñ 0 and K Z Z =    k(z 1 , z1 ) • • • k(z 1 , zñ 0 ) . . . . . . . . . k(z ñ0 , z1 ) • • • k(z ñ0 , zñ 0 )    ∈ R ñ0 ×ñ 0 .Set P Z h * = ñ0 j=1 a j k zj . Then we can see that(I -ηS Z SZ ) ñ0 j=1 ajkz j = ñ0 r=1 ar - η n n 0 i=1 ñ0 j=1 ajk(zj, zi)bir kz r = ñ0 r=1 ar - η n [B BK Z Z a]r kz r = ñ0 r=1 I - η n B BK Z Z a r kz r where [•] r is the rth component of the given vector and a = [a 1 , • • • , a ñ0 ] . Note that K Z Z is invertible since v K Z Z v = 0 implies v = 0.We can also see thatK ZZ = BK Z Z B where KZZ =    k(z 1 , z 1 ) • • • k(z 1 , z n 0 ) . . . . . . . . . k(z n 0 , z 1 ) • • • k(z n 0 , z n 0 )    ∈ R n 0 ×n 0 . So K 1/2 Z Z B BK 1/2 Z Z = BK Z Z B ≤ BK Z Z B F ≤ nκ 2 .Thus, 0< η n K 1/2 Z Z B BK 1/2 Z Z < I and (I -ηS Z SZ ) t PZ h * = ñ0 r=1 K -1/2 Z Z I - η n K 1/2 Z Z B BK 1/2 Z Z t K 1/2 Z Z a r kz r → 0 as t → ∞. Similarly, we get (I -ηS Z S Z ) t u 1 → 0 as t → ∞. Therefore, we attain G t u 1 → P Z h * as t → ∞ for any u 1 ∈ W .From this fact, DCL-KR has the recurrence relationEQUATIONwheref t = f i,t for any i = 1, • • • , m and T k,Xi := I -ηT k,Xi for i = 1, • • • , m.Then we obtain a closed formft = PZ m i=1 ni n T E k,X i t f0 + t-1 j=0 PZ m i=1 ni n T E k,X i j PZ m i=1 ni n η E-1 s=0 T s k,X i S D i yi.We first computef * 0 -   PZ m i=1 ni n T E k,X i t f * 0 + t-1 j=0 PZ m i=1 ni n T E k,X i j PZ m i=1 ni n η E-1 s=0 T s k,X i T k,X i f * 0   . From η E-1 s=0 T s k,X i T k,X i = η E-1 s=0 (I -ηT k,X i ) s T k,X i = I -(I -ηT k,X i ) E ,we havef * 0 -   PZ m i=1 ni n T E k,X i t f * 0 + t-1 j=0 PZ m i=1 ni n T E k,X i j PZ m i=1 ni n η E-1 s=0 T s k,X i T k,X i f * 0   = f * 0 -   PZ m i=1 ni n T E k,X i t f * 0 + t-1 j=0 PZ m i=1 ni n T E k,X i j PZ I - m i=1 ni n T E k,X i f * 0   = I + m i=1 ni n PZ T E k,X i + • • • + m i=1 ni n PZ T E k,X i t-1 (I -PZ )f * 0 = (I -PZ )f * 0 + I + • • • + m i=1 ni n PZ T E k,X i t-2 PZ m i=1 ni n T E k,X i -I (I -PZ )f * 0 .where the last equality follows from P Z (I -P Z ) = 0. Thus, we obtain the equalityEQUATIONA.2.2 Norm Bound of First Term in (7) We first bound the norm of the first term in (7) asEQUATIONwhere λ > 0. The first term in ( 8) is bounded asT 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i t PZ f * 0 H k ≤ T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i t PZ (T k,X + λI) r-1/2 • (T k,X + λI) -(r-1/2) T r-1/2 k,ρx g * 0 H k .Note that(T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i t PZ (T k,X + λI) r-1/2 ≤ (T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i PZ t/2r 2r by Lemma A.8. Set A i = T E k,Xi ⇔ T k,Xi = 1 η (I -A 1/E i). We observe that(T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i PZ t/2r 2r = PZ m i=1 ni n T E k,X i PZ t/2r PZ (T k,X + λI)PZ PZ m i=1 ni n T E k,X i PZ t/2r r ≤   1 η m i=1 ni n PZ AiPZ t/2r I - n i=1 ni n PZ AiPZ m i=1 ni n PZ AiPZ t/2r + λ   rwhere the equality follows from 0≤ A i ≤ I ⇒ A 1/E i ≥ A i and I ≥ P Z . Since sup x∈[0,1] x t/r (1- x) = r t+r • ( t t+r ) t/r, we attain the inequalityEQUATIONNext, Lemma A.8 gives(T k,X + λI) -(r-1/2) T r-1/2 k,ρx ≤ (T k,X + λI) -(r-1/2) (T k,ρx + λI) r-1/2 ≤ (T k,X + λI) -1 (T k,ρx + λI) r-1/2 and T 1/2 k,ρx (T k,X + λI) -1/2 ≤ (T k,ρx + λI) 1/2 (T k,X + λI) -1/2 ≤ (T k,ρx + λI)(T k,X + λI) -1 1/2 .By Lemma A.10,EQUATIONholds with confidence at least 1 -δ where δ ∈ (0, 1). Combining ( 9) and ( 10) and applying r t+r ≤ 1 t and ( t t+r) t/r ≤ 1 2 yield T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i t PZ f * 0 H k ≤ r t + r • 1 η t t + r t/r + λ r g * 0 H k 2 + 2 2κ 2 nλ + 4κ 2 N (λ) nλ log(2/δ) 2 r ≤ R 1 2ηt + λ r 2 + 2 2κ 2 nλ + 4κ 2 N (λ) nλ 2 r (log(4/δ)) 2rwith confidence at least 1 -δ where δ ∈ (0, 1). Therefore, putting λ = n -1 2r+s yieldsE   T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i t PZ f * 0 H k   ≤ 1 2ηt + n -1 2r+s r R • 4Γ(2r + 1) 2 + 2(2κ 2 + 2κ C s ) 2 r 1 t + n -1 2r+s r .Here, we apply the fact that EA = ∞ 0 P(A ≥ t) dt for A ≥ 0. We next turn to bound the second term in (8) . Note that(T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i t = m i=1 ni n T E k,X i PZ t (T k,X + λI) PZ m i=1 ni n T E k,X i t 1/2 ≤ m i=1 ni n T E k,X i • PZ m i=1 ni n T E k,X i PZ t-1 PZ (T k,X + λI)PZ PZ m i=1 ni n T E k,X i PZ t-1 1/2 ≤ PZ m i=1 ni n T E k,X i PZ t-1 PZ (T k,X + λI)PZ PZ m i=1 ni n T E k,X i PZ t-1 1/2 . Set A i = T E k,Xi. Using a similar argument as before givesPZ m i=1 ni n T E k,X i PZ t-1 PZ (T k,X + λI)PZ PZ m i=1 ni n T E k,X i PZ t-1 1/2 ≤ 1 η(2t -1) + λ 1/2 ≤ 1 ηt + λ 1/2 .Since Z and X are independent, we haveE T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i t (I -PZ )f * 0 H k ≤ 1 ηt + λ 1/2 E T 1/2 k,ρx (T k,X + λI) -1/2 • E (I -PZ )f * 0 H k .We already see thatT 1/2 k,ρx (T k,X + λI) -1/2 ≤ 2 + 2 2κ 2 nλ + 4κ 2 N (λ) nλ log(2/δ) 2 1/2 ≤ 2 + 2 2κ 2 nλ + 4κ 2 N (λ) nλ 2 1/2 log(4/δ)holds with confidence at least 1 -δ where δ ∈ (0, 1) and soE T 1/2 k,ρx (T k,X + λI) -1/2 ≤ 4 2 + 2(2κ 2 + 2κ C s ) 2 1/2by putting λ = n -1 2r+s as before. The remaining part is to boundE (I -P Z )f * 0 H k . Applying Lemma A.9 yields (I -P Z )f * 0 H k ≤ λ 1/2 0 (T k,Z + λ 0 I) -1/2 T r-1/2 k,ρx g * 0 H k where λ 0 > 0. Then Lemma A.8 gives λ 1/2 0 (T k,Z + λ0I) -1/2 T r-1/2 k,ρx g * 0 H k ≤ Rλ 1/2 0 (T k,Z + λ0I) -(1-r) (T k,Z + λ0I) -(r-1/2) T r-1/2 k,ρx ≤ Rλ r-1/2 0 (T k,Z + λ0I) -1/2 T 1/2 k,ρx 2r-1 .From dρx d ρx ≤ B, we obtainT 1/2 k,ρx (T k,Z + λ0I) -1/2 = ιρ x (T k,Z + λ0I) -1/2 ≤ B 1/2 ιρ x (T k,Z + λ0I) -1/2 = B 1/2 T 1/2 k, ρx (T k,Z + λ0I) -1/2 .Set λ 0 = 128(κ 2 + 1) 2 (log n 0 ) 3 /n 0 where we assume n is sufficiently large such that λ 0 ≤ 1 and N ρx (λ 0 ) ≥ 1 for n 0 ≥ n 1 2r+s (log n) 3 . By Lemma A.8 and Lemma A.12,T 1/2 k, ρx (T k,Z + λ0I) -1/2 ≤ (T k, ρx + λ0I) 1/2 (T k,Z + λ0I) -1/2 ≤ √ 2holds with confidence at least 1 -δ where δ ∈ [4 exp(-1/4(κ 2 + 1)B 0 ), 1) andB0 = 1 + log Nρ x (λ0) λ0n0 + 1 + log Nρ x (λ0) λ0n0 .Since(I -PZ )f * 0 H k ≤ f * 0 H k = T r-1/2 k,ρx g * 0 H k ≤ T k,ρx r-1/2 g * 0 H k ≤ Rκ 2r-1 ,we haveE (I -PZ )f * 0 H k ≤ Rλ r-1/2 0 B r-1/2 2 r-1/2 + Rκ 2r-1 • 4 exp - 1 4(κ 2 + 1)B0.From n 0 ≥ κ 2 e, we getB0 ≤ log κ 2 e + log n0 128(κ 2 + 1) 2 (log n0) 3 + log κ 2 e + log n0 128(κ 2 + 1) 2 (log n0) 3 ≤ 2 log n0 128(κ 2 + 1) 2 (log n0) 3 + 2 log n0 128(κ 2 + 1) 2 (log n0) 3 ≤ 1 4(κ 2 + 1) log n0 and so Rκ 2r-1 • 4 exp - 1 4(κ 2 +1)B0 ≤ 4Rκ 2r-1 • 1 n0 . Therefore, E T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 PZ m i=1 ni n T E k,X i t (I -PZ )f * 0 H k ≤ 1 ηt + n -1 2r+s 1/2 4 2 + 2(2κ 2 + 2κ C s ) 2 1/2 Rλ r-1/2 0 B r-1/2 2 r-1/2 + 4Rκ 2r-1 • 1 n0 B r-1/2 1 t + n -1 2r+s 1/2 n -r-1/2 2r+swhere the last inequality comes from n 0 ≥ n 1 2r+s (log n) 3 .A.2.3 Norm Bound of Second Term in ( 7)Set P =      E-1 s=0 (I -ηSD 1 S D 1 ) s 0 • • • 0 0 E-1 s=0 (I -ηSD 2 S D 2 ) s • • • 0 . . . . . . . . . . . . 0 0 • • • E-1 s=0 (I -ηSD m S Dm ) s      .Note thatEQUATIONThen the second term in (7) becomesιρ x t-1 j=0 PZ m i=1 ni n T E k,X i j PZ m i=1 ni n η E-1 s=0 T s k,X i S D i (yi -SD i f * 0 ) = ιρ x t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ).We can see thatιρ x t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) L 2 ρx ≤ T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) H k ≤ T 1/2 k,ρx (T k,X + λI) -1/2   T 1/2 k,X t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) H k +λ 1/2 t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) H k   .We first bound the expectation of the first term in the above. By the Cauchy-Schwartz inequality, we haveE   T 1/2 k,ρx (T k,X + λI) -1/2 T 1/2 k,X t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) H k   ≤ E T 1/2 k,ρx (T k,X + λI) -1/2 2 1/2   E T 1/2 k,X t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) 2 H k   1/2 .Observe thatT 1/2 k,X t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) H k = SD t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) 2 and SD t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) = P -1/2 (I -(I -ηP 1/2 SDPZ S D P 1/2 ) t )P 1/2 (y -SDf * 0 ) = (I -(I -ηSDPZ S D P ) t )(y -SDf * 0 ). Using E(y -S D f * 0 )(y -S D f * 0 ) ≤ γ 2 I,we haveE (I -(I -ηSDPZ S D P ) t )(y -SDf * 0 ) 2 2 |X, Z = 1 n tr (I -(I -ηSDPZ S D P ) t )E (y -SDf * 0 )(y -SDf * 0 ) (I -(I -ηSDPZ S D P ) t ) ≤ γ 2 n (I -(I -ηSDPZ S D P ) t ) 2 HS ≤ γ 2 E n (I -(I -ηP 1/2 SDPZ S D P 1/2 ) t ) 2 HSwhere the last inequality follows from the fact that AB HS ≤ A B HS , AB HS ≤ A HS B , andEQUATIONwhich follows from (11) , we can see that0 ≤ λi(ηP 1/2 SDPZ S D P 1/2 ) ≤ λi(ηP 1/2 SDS D P 1/2 ) ≤ 1 ⇒ 0 ≤ λi(I -(I -ηP 1/2 SDPZ S D P 1/2 ) t ) ≤ λi(I -(I -ηP 1/2 SDS D P 1/2 ) t ) ≤ 1where λ i (•) is the ith largest eigenvalue of a given operator. Therefore,γ 2 E n (I -(I -ηP 1/2 SDPZ S D P 1/2 ) t ) 2 HS ≤ γ 2 E n (I -(I -ηP 1/2 SDS D P 1/2 ) t ) 2 HS .Using (12) and1 ∧ u 2 ≤ 1 ∧ u for u ≥ 0 lead to λi(I -(I -ηP 1/2 SDS D P 1/2 ) t ) 2 = (1 -(1 -ηλi(P 1/2 SDS D P 1/2 )) t ) 2 ≤ 1 ∧ (η 2 t 2 λi(P 1/2 SDS D P 1/2 ) 2 ) ≤ 1 ∧ (ηtλi(P 1/2 SDS D P 1/2 )) ≤ 1 ∧ (ηtE λi)where λ1 ≥ • • • ≥ λn are eigenvalues of S D S D , the first inequality comes from the Bernoulli inequality, and the last inequality follows from the fact that P ≤ E. We defineR( ) = 1 n n i=1 ( λi ∧ 2 ).Thenγ 2 E n (I -(I -ηP 1/2 SDS D P 1/2 ) t ) 2 HS ≤ γ 2 ηtE 2 • R 1 √ ηtE 2 .Similarly as in Appendix A.2.2, putting λ = n -1 2r+s givesE T 1/2 k,ρx (T k,X + λI) -1/2 2 ≤ 2 + 4Γ(3)(2κ 2 + 2κ C s ) 2 .Therefore, the Cauchy-Schwartz inequality gives a bound asE   T 1/2 k,ρx (T k,X + λI) -1/2 T 1/2 k,X t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) H k   ≤ (2 + 4Γ(3)(2κ 2 + 2κ C s ) 2 )γ 2 ηtE 2 • ER 1 √ ηtE 2 1/2.We now bound the expectation ofλ 1/2 T 1/2 k,ρx (T k,X + λI) -1/2 t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) H k .By the Cauchy-Schwartz inequality and the same argument as before, we haveE   λ 1/2 T 1/2 k,ρx (T k,X + λI) -1/2 t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) H k   ≤ (2 + 4Γ(3)(2κ 2 + 2κ C s ) 2 ) 1/2   λ • E t-1 j=0 PZ -ηPZ S D P SDPZ j ηPZ S D P (y -SDf * 0 ) 2 H k   1/2 .Also, the same argument as before yields To bound E[tr(A)], note thatE t-1 j=0 PZ -ηPZ S D P SDPZ j ηPZ S D P (y -SDf * 0 ) 2 H k = 1 n E (y -SDf * 0 ) A(y -SDf * 0 ) ≤ γ 2 n E[tr(A)]tr(A) ≤ E • tr   ηP 1/2 SDPZ t-1 j=0 I -ηPZ S D P SDPZ j 2 ηPZ S D P 1/2   = ηE • tr   ηP 1/2 SDPZ S D P 1/2 t-1 j=0 (I -ηP 1/2 SDPZ S D P 1/2 ) j 2   . Let B = ηP 1/2 S D P Z S D P 1/2 . Then 0 ≤ B ≤ I and ηE • tr   ηP 1/2 SDPZ S D P 1/2 t-1 j=0 (I -ηP 1/2 SDPZ S D P 1/2 ) j 2   = ηE n i=1 λi(B) t-1 j=0 (1 -λi(B)) j 2 = ηE n i=1 1 λi(B) (1 -(1 -λi(B)) t ) 2 ≤ ηE n i=1 1 λi(B) ∧ (t 2 λi(B)) ≤ ηE n i=1 t ∧ (t 2 λi(B))where the first inequality follows from 1 -x t ≤ 1 ∧ t(1 -x) and the second inequality follows from 1/x ∧ t 2 x ≤ t ∧ t 2 x for all t ≥ 0 and x ∈ [0, 1]. From the fact thatλi(B) ≤ η P 1/2 2 λi(SDPZ S D ) ≤ ηEλi(SDS D ) = ηE λi,we haveηE n i=1 t ∧ (t 2 λi(B)) ≤ ηE n i=1 t ∧ (ηt 2 E λi) = nη 2 t 2 E 2 • R 1 √ ηtE 2 .Therefore, we obtainE   λ 1/2 T 1/2 k,ρx (T k,X + λI) -1/2 t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) H k   ≤ (2 + 4Γ(3)(2κ 2 + 2κ C s ) 2 )(λγ 2 η 2 t 2 E 2 ) • ER 1 √ ηtE 2 1/2 .In conclusion, we have an upper bound of the norm of the second term in (7) asE ιρ x t-1 j=0 PZ -ηPZ S D P SD j ηPZ S D P (y -SDf * 0 ) L 2 ρx ≤ (2 + 4Γ(3)(2κ 2 + 2κ C s ) 2 )γ 2 ηtE 2 • ER 1 √ ηtE 2 1/2 + (2 + 4Γ(3)(2κ 2 + 2κ C s ) 2 )(λγ 2 η 2 t 2 E 2 ) • ER 1 √ ηtE 2 1/2 t 1/2 + n -1/2 2r+s t • ER 1 √ ηtE 2 1/2by taking λ = n -1 2r+s . We will bound ER1 √ ηtE 2 in Appendix A.2.5.A.2.4 Norm Bound of Third and Last Term in (7) Note that 11) and 0 ≤ P ≤ EI, we have 0 ≤ I -ιρ x I + • • • + m i=1 ni n PZ T E k,X i t-2 PZ I - m i=1 ni n T E k,X i (I -PZ )f * 0 L 2 ρx ≤ T 1/2 k,ρx (T k,X + λI) -1/2 • (T k,X + λI) 1/2 I + • • • + m i=1 ni n PZ T E k,X i t-2 PZ I - m i=1 ni n T E k,X i 1/2 • I - m i=1 ni n T E k,X i 1/2 (I -PZ )f * 0 H k where 0 < λ ≤ 1. From (m i=1 ni n T E k,Xi = ηS D P S D ≤ ηES D S D = ηET k,X .Using this fact, we find that(T k,X + λI) 1/2 I + • • • + m i=1 ni n PZ T E k,X i t-2 PZ I - m i=1 ni n T E k,X i 1/2 ≤ (ηE) 1/2 T 1/2 k,X I + • • • + m i=1 ni n PZ T E k,X i t-2 PZ T 1/2 k,X +λ 1/2 I + • • • + m i=1 ni n PZ T E k,X i t-2 PZ T 1/2 k,X .To bound this, we first observe thatT 1/2 k,X I + • • • + m i=1 ni n PZ T E k,X i t-2 PZ T 1/2 k,X ≤ t-2 j=0 T 1/2 k,X PZ m i=1 ni n PZ T E k,X i PZ j PZ T 1/2 k,X ≤ t-2 j=0 m i=1 ni n PZ T E k,X i PZ j/2 PZ T k,X PZ m i=1 ni n PZ T E k,X i PZ j/2 ≤ 1 η 1 + t-2 j=1 1 jwhere the last inequality follows by a similar calculation as before. On the other hand, we haveI + • • • + m i=1 ni n PZ T E k,X i t-2 PZ T 1/2 k,X ≤ t-2 j=0 m i=1 ni n PZ T E k,X i PZ j PZ T 1/2 k,X = t-2 j=0 m i=1 ni n PZ T E k,X i PZ j PZ T k,X PZ m i=1 ni n PZ T E k,X i PZ j 1/2 ≤ 1 √ η 1 + t-2 j=1 1 √ 2jby the same argument. Using a simple calculation, we get1 η 1 + t-2 j=1 1 j ≤ 1 ηE (2 + log t)E and 1 √ η 1 + t-2 j=1 1 √ 2j ≤ 1 √ η + 1 √ 2η (2 √ t -2 -1) ≤ 1 (ηE) 1/2 • √ 6tE.Note that the norm of the third term in ( 7) is bounded asιρ x (I -PZ )f * 0 L 2 ρx ≤ T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 (I -PZ )f * 0 H k .Therefore, the norm of the sum of the third and last terms in ( 7) is bounded by(1 + 2E + E log t + 6ηtλE) T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 (I -PZ )f * 0 H k . To bound (T k,X + λI) 1/2 (I -P Z )f * 0 H k , observe that (T k,X + λI) 1/2 (I -PZ )f * 0 H k ≤ (T k,X + λI) 1/2 (T k,ρx + λI) -1/2 (T k,ρx + λI) 1/2 (I -PZ )f * 0 H k and (T k,ρx + λI) 1/2 (I -PZ )f * 0 2 H k = ιρ x (I -PZ )f * 0 2 L 2 ρx + λ (I -PZ )f * 0 2 H k ≤ B ιρ x (I -PZ )f * 0 2 L 2 ρx + λ (I -PZ )f * 0 2 H k = B (T k, ρx + λI) 1/2 (I -PZ )f * 0 2 H k .Under Assumption 3.3, we haveB 1/2 (T k, ρx + λI) 1/2 (I -PZ )f * 0 H k ≤ B 1/2 (T k, ρx + λI) 1/2 (I -PZ )(T k, ρx + λI) r-1/2 (T k, ρx + λI) -(r-1/2) T r-1/2 k,ρx g * 0 H k . Since (T k, ρx + λI) -(r-1/2) T r-1/2 k,ρx ≤ (T k, ρx + λI) -1/2 T 1/2 k,ρx 2r-1 = T 1/2 k,ρx (T k, ρx + λI) -1/2 2r-1which follows from Lemma A.8 andT 1/2 k,ρx (T k, ρx + λI) -1/2 = ιρ x (T k, ρx + λI) -1/2 ≤ B 1/2 ιρ x (T k, ρx + λI) -1/2 ≤ B 1/2 T 1/2 k, ρx (T k, ρx + λI) -1/2 ≤ B 1/2 ,we have (T k, ρx + λI) -(r-1/2) T r-1/2 k,ρx ≤ B r-1/2 . On the other hand,(T k, ρx + λI) 1/2 (I -PZ )(T k, ρx + λI) r-1/2 ≤ (T k, ρx + λI) 1/2 (I -PZ ) (I -PZ ) 2r-1 (T k, ρx + λI) r-1/2 ≤ (T k, ρx + λI) 1/2 (I -PZ ) 2rby Lemma A.8. Therefore,T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 (I -PZ )f * 0 H k ≤ RB r (T k,ρx + λI) 1/2 (T k,X + λI) -1/2 • (T k,X + λI) 1/2 (T k,ρx + λI) -1/2 (T k, ρx + λI) 1/2 (I -PZ ) 2r .Since X and Z are independent, we haveE T 1/2 k,ρx (T k,X + λI) -1/2 (T k,X + λI) 1/2 (I -PZ )f * 0 H k ≤ RB r E (T k,ρx + λI) 1/2 (T k,X + λI) -1/2 (T k,X + λI) 1/2 (T k,ρx + λI) -1/2 • E (T k, ρx + λI) 1/2 (I -PZ ) 2r .By Lemma A.8 and Lemma A.10,(T k,ρx + λI) 1/2 (T k,X + λI) -1/2 ≤ 2 + 2 2κ 2 nλ + 4κ 2 N (λ) nλ log(2/δ) 2 1/2holds with confidence at least 1 -δ where δ ∈ (0, 1). Also, by Lemma A.8 and Lemma A.11(T k,X + λI) 1/2 (T k,ρx + λI) -1/2 ≤ 1 + 2κ 2 nλ + 4κ 2 N (λ) nλ log(2/δ) 1/2holds with confidence at least 1 -δ where δ ∈ (0, 1). Thus,(T k,ρx + λI) 1/2 (T k,X + λI) -1/2 (T k,X + λI) 1/2 (T k,ρx + λI) -1/2 ≤ 2 + 2 2κ 2 nλ + 4κ 2 N (λ) nλ 2 1/2 1 + 2κ 2 nλ + 4κ 2 N (λ) nλ 1/2(log(4/δ)) 3/2 with confidence at least 1 -δ where δ ∈ (0, 1). Set λ = 128(κ 2 + 1) 2 n -1 2r+s where n is sufficiently large such that λ ≤ 1 and N ρx (λ) ≥ 1. ThenE (T k,ρx + λI) 1/2 (T k,X + λI) -1/2 (T k,X + λI) 1/2 (T k,ρx + λI) -1/2 ≤ 4Γ (2.5) (2 + 2(2κ 2 + 2κ C s ) 2 ) 1/2 (1 + 2κ 2 + 2κ C s ) 1/2 1.We now bound E (T k, ρx + λI) 1/2 (I -P Z ) 2r . By Lemma A.9, we have(T k, ρx + λI) 1/2 (I -PZ ) 2r ≤ λ r (T k, ρx + λI) 1/2 (T k,Z + λI) -1/2 2r .By Lemma A.12, (T k, ρx + λI) 1/2 (T k,Z + λI) -1/2 ≤ √ 2 with confidence at least 1 -4 exp(-1/4(κ 2 + 1)B 0 ) whereB0 = 1 + log Nρ x (λ) λn0 + 1 + log Nρ x (λ) λn0 .Also, (T k, ρx + λI) 1/2 (I -P Z ) ≤ (κ 2 + 1) 1/2 almost surely. Thus,E (T k, ρx + λI) 1/2 (I -PZ ) 2r ≤ 2 r λ r + (κ 2 + 1) r • 4 exp - 1 4(κ 2 + 1)B0 .Note thatB0 ≤ log κ 2 e + log(1/λ) λn0 + log κ 2 e + log(1/λ) λn0 ≤ 1 4(κ 2 + 1) log n and so (κ 2 + 1) r • 4 exp - 1 4(κ 2 +1)B0 ≤ 4(κ 2 + 1) r • 1 n . Therefore, E (T k, ρx + λI) 1/2 (I -PZ ) 2r ≤ 2 r 128 r (κ 2 + 1) 2r • n -r 2r+s + 4(κ 2 + 1) r • n -1 n -r 2r+s .We can conclude thatE -ιρ x (I -PZ )f * 0 + ιρ x I + • • • + m i=1 ni n PZ T E k,X i t-2 PZ I - m i=1 ni n T E k,X i (I -PZ )f * 0 L 2 ρx ≤ (1 + 2E + E log t + 6ηtλE)RB r • 4Γ (2.5) (2 + 2(2κ 2 + 2κ C s ) 2 ) 1/2 (1 + 2κ 2 + 2κ C s ) 1/2 • 2 r 128 r (κ 2 + 1) 2r • n -r 2r+s + 4(κ 2 + 1) r • n -1 B r (1 + log t + t 1/2 n -1/2 2r+s )n -r 2r+s .For convenience, we abuse the notationD = {(x 1 , y 1 ), • • • , (x n , y n )} and X = {x 1 , • • • , x n }.Define the local empirical Rademacher complexityQn( ) = E   sup g H k ≤1, g L 2 ρx,n ≤ 1 n n i=1 wig(x i ) X  and the local population Rademacher complexityQ n ( ) = E   sup g H k ≤1, g L 2 ρx ≤ 1 n n i=1 wig(x i )  where w 1 , • • • , w n are independent Rademacher random variables and ρ x,n = 1 n n i=1 δ x i . We also defineR( ) = 1 n ∞ i=1 λi ∧ 2where λ 1 ≥ λ 2 ≥ • • • ≥ 0 are eigenvalues of T k,ρx . We recall the following well-known property. Lemma A.1 ([46], [60] ). We haveQ n ( ) ≤ √ 2 • R( ) for > 0.We can prove the following lemma using a similar argument as in [46] . Lemma A.2. There is an absolute constant c > 0 which satisfies that for every > 0,c • R( ) ≤ Qn( ).Proof of Lemma A.2. We divide the proof into three parts.Part 1. Since T k,X = S D S D , λ1 ≥ λ2 ≥ • • • ≥ λn ≥ 0 are eigenvalues of T k,X .For convenience, set λi = 0 for i > n and define n ≤ n such that λn > 0 and λn+1 = 0. Choose an orthonormal basis{ ψi } ∞ i=1 of H k such that ψi is an eigenvector of T k,X corresponding to λi . Then ψi , ψj L 2 ρx,n = S D ψi , S D ψj 2 = T k,X ψi , ψj H k = δ {i=j} λi . We will show that kx = n i=1 ψi(x) ψi where x ∈ {x 1 , • • • , x n }. Let W 1 be the subspace of H k spanned by { ψi : i = 1, • • • , n} and W 2 be the subspace of H k spanned by {k x i : i = 1, • • • , n}. Observe that W ⊥ 1 = ker T k,X and W ⊥2 ⊂ ker T k,X by the reproducing property. Thus, W 1 ⊂ W 2 . Conversely, choose a basis{k xi : i = 1, • • • , ñ } ⊂ {k x i : i = 1, • • • , n} of W 2 .Then, using a similar argument as in Appendix A.2.1 implies that there exists a matrixB =    b11 • • • b 1ñ . . . . . . . . . bn1 • • • b nñ    ∈ R n×ñ such that k x i = ñ j=1 b ij k xj . Then K X X = BK X X where K X X =     k(x 1 , x1 ) • • • k(x 1 , xñ ) . . . . . . . . . k(x n , x1 ) • • • k(x n , xñ )     ∈ R n×ñ and K X X =     k(x 1 , x1 ) • • • k(x 1 , xñ ) . . . . . . . . . k(x ñ , x1 ) • • • k(x ñ , xñ )     ∈ R ñ ×ñ .Since K X, X and B B are invertible,T k,X   ñ i=1 [nK -1 X, X (B B) -1 b]ik xi   = ñ i=1 bik xi for any b = [b 1 , • • • , b ñ ] ∈ R ñ where [•] r is the rth component of the given vector. Therefore, W 2 ⊂ ran T k,X = (ker T k,X ) ⊥ = W 1 and so W 1 = W 2 .From this fact, we can see thatk x i = n r=1 a r ψr for some a 1 , • • • , a n ∈ R. Then a r = k x i , ψr H k = ψr (x i ) for all r = 1, • • • , n and so we are done. Note that k x = ∞ i=1 ψi (x) ψi where x ∈ {x 1 , • • • , x n } since ψi (x) = 0 for x ∈ {x 1 , • • • , x n } and i > n. Part 2. Define F := h : h H k ≤ 1 and h L 2 ρx,n ≤ = ∞ i=1 hi ψi : ∞ i=1 h 2 i ≤ 1 and n i=1 λih 2 i ≤ 2 and E := ∞ i=1 hi ψi : ∞ i=1 λi λi ∧ 2 h 2 i ≤ 1 where 0 0 = 1. Then E ⊂ F since ∞ i=1 h 2 i ∨ ∞ i=1 λi 2 h 2 i ≤ ∞ i=1 1 ∨ λi 2 h 2 i = ∞ i=1 λi λi ∧ 2 h 2 i ≤ 1 for h = ∞ i=1 h i ψi ∈ E. Thus, E sup h∈E n i=1 ih(x i ) 2 X ≤ E sup h∈F n i=1 ih(x i ) 2 Xwhere 1 , • • • , n are i.i.d. Rademacher variables. By the reproducing property,n i=1 ih(x i ) = h, n i=1 ikxi H k = h, n i=1 i n j=1 ψj(x i ) ψj H k = n j=1hjn i=1 i ψj(x i ) = ∞ j=1 λj λj ∧ 2 hj ψj, n j=1 λj ∧ 2 λj n i=1 i ψj(x i ) ψj H k where h = ∞ i=1 h i ψi . Thus, sup h∈E n i=1 ih(x i ) 2 = n j=1 λj ∧ 2 λj n i=1 i ψj(x i ) ψj 2 H k . Since n j=1 λj ∧ 2 λj n i=1 i ψj(x i ) ψj 2 H k = n j=1 λj ∧ 2 λj n i=1 i ψj(x i ) 2 ,we haveE sup h∈E n i=1 ih(x i ) 2 X = E n j=1 λj ∧ 2 λj n i=1 i ψj(x i ) 2 X = n j=1 λj ∧ 2 λj n i=1 ψj(x i ) 2 = n n j=1 λj ∧ 2 .Therefore,n n j=1 λj ∧ 2 ≤ E sup h∈F n i=1 ih(x i ) 2 X 1/2 .Part 3. By Khintchine's inequality,E sup h∈F n i=1 ih(x i ) X ≥ sup h∈F E n i=1 ih(x i ) X ≥ 1 √ 2 sup h∈F n i=1 h(x i ) 2 1/2 = n 2 . Set Z = g( 1 , • • • , n ) where g(t1, • • • , tn) = sup h∈F n i=1 tih(x i ) .By Remark A.14, g is convex and sup h∈F (n i=1 h(x i ) 2 ) 1/2 = √ n -Lipschitz on [-1, 1] n. By Lemma A.13, we haveP (Z -E [Z|X] ≥ tE [Z|X] |X) ≤ exp - t 2 16 2 n E [Z|X] 2 ≤ exp - t 2 32 .FromE Z 2 |X = E Z 2 1 {Z<E[Z|X]} |X + ∞ m=0 E Z 2 1 {(m+1)E[Z|X]≤Z<(m+2)E[Z|X]} |X ≤ E[Z|X] 2 + ∞ m=0 (m + 2) 2 E[Z|X] 2 P(Z ≥ (m + 1)E[Z|X]|X) ≤ E[Z|X] 2 1 + ∞ m=0 (m + 2) 2 exp - m 2 32 ,we haveE[Z|X] ≥ c • E Z 2 |X 1/2 ≥ c • n n j=1 λj ∧ 2 where c = 1 + ∞ m=0 (m + 2) 2 exp -m 2 32 -1/2is an absolute constant. Therefore,c • 1 n n j=1 λj ∧ 2 ≤ E sup h∈F 1 n n i=1 ih(x i ) X .We set the population radius asn = inf ≥ 0 : Q n ( ) ≤ 1+2r 16κ .We also define˜ n = inf ≥ 0 : R( ) ≤ 1+2r 16 √ 2κ .By Lemma A.1, we have n ≤ ˜ n . We can easily see that R, R, Q n , and Q n are increasing functions.The following lemma can be shown by a similar argument as in [2] .Lemma A.3. If g : [0, ∞) → [0, ∞)is a function such that g is non-decreasing and r → g(r)/r is non-increasing, then g is continuous on (0, ∞).R( ) = 1 n n i=1 1 ∧ λi 2 and R( ) = 1 n ∞ i=1 1 ∧ λi 2 ,→ R( )/ and → R( )/ are non-increasing and so R and R are continuous.Lemma A.4. → Q n ( )/ and → Q n ( )/ are non-increasing. In particular, Q n is continuous and n < ∞.Proof. From the fact thatQn( ) = 1 E   sup g H k ≤1, g L 2 ρx,n ≤ 1 n n i=1 wig(xi) X   = E   sup g H k ≤1/ , g L 2 ρx,n ≤1 1 n n i=1 wig(xi) X   ,we can easily see that → Q n ( )/ is non-increasing. Similarly, we can show that → Q n ( )/ is non-increasing. Note that lim →0 + Qn( ) > 0 and lim→0 + Q n ( ) > 0.Also, we can observe thatlim →∞ Qn( ) = 0 and lim →∞ Q n ( ) = 0.Since → 2r /16κ is increasing, goes 0 as → 0 + , and goes ∞ as → ∞, we can conclude that n < ∞.Similarly, we have ˜ n < ∞. In fact, we can find the lower and the upper bound of ˜ n under Assumption 3.2. Lemma A.5. We have2 9/(4r+2s) κ 1/(2r+s) c s/(4r+2s) s ∧ c 1/2 s s s + 2 1/2s n -1 4r+2s ≤ ˜ n ≤ 2 9/(4r+2s) κ 1/(2r+s) 2 -s 1 -s 1/(4r+2s) C s/(4r+2s) s n -1 4r+2s . Proof. Since c s i -1/s ≤ λ i ≤ C s i -1/s , we have 1 n ∞ j=1 (csj -1/s ) ∧ 2 ≤ R( ) ≤ 1 n ∞ j=1 (Csj -1/s ) ∧ 2 .We first consider the lower bound of ˜ n . We first observe that1 n ∞ j=1 (csj -1/s ) ∧ 2 = 1 n    cs 2 s 2 + ∞ j= cs 2 s +1 csj -1/s   . Set = 2 9/(4r+2s) κ 1/(2r+s) c s/(4r+2s) s ∧ c 1/2 s s s + 2 1/2s n -1 4r+2s .Note thatcs cs 2 s -1/s ≥ 2 and s 1 -s cs 2 s + 1 (s-1)/s ≥ cs 2 s -1/shold. The first formula is trivial. To show the second formula, we observe that the functionu(t) = ( s 1-s ) s • t (1+t) 1-s is an increasing function. Thus, for t ≥ 2/s u(t) ≥ u 2 s = 2 (1 -s) s (s + 2) 1-s ≥ 2 2 -2s 2 ≥ 1.Here, we apply an elementary inequality:a s b 1-s ≤ sa + (1 -s)b ∀a, b > 0. Since ≤ c 1/2 s s s + 2 1/2s ⇒ cs 2 s ≥ cs 2 s -1 ≥ 2 s , putting t = cs 2 s gives s 1 -s s • cs 2 s (1 + cs 2 s ) 1-s ≥ 1and so the second formula holds. Therefore,∞ j= cs 2 s +1 csj -1/s ≥ cs ∞ cs 2 s +1 1 t 1/s dt = scs 1 -s cs 2 s + 1 (s-1)/s ≥ cs cs 2 s -1/s ≥ 2holds and so we1 n    cs 2 s 2 + ∞ j= cs 2 s +1 csj -1/s    ≥ 1 n cs 2 s + 1 2 ≥ c s/2 s n 1/2 1-s ≥ 1+2r 16 √ 2κwhere the last inequality follows from ≤ 2 9/(4r+2s) κ 1/(2r+s) c s/(4r+2s) s n -1 4r+2s . We can conclude that˜ n ≥ 2 9/(4r+2s) κ 1/(2r+s) c s/(4r+2s) s ∧ c 1/2 s s s + 2 1/2s n -1 4r+2sby Lemma A.4. We now derive the upper bound of ˜ n . Note that1 n ∞ j=1 (Csj -1/s ) ∧ 2 = 1 n    Cs 2 s 2 + ∞ j= Cs 2 s +1 Csj -1/s   . Since ∞ j= Cs 2 s +1 Csj -1/s - ∞ Cs 2 s +1 Cs t 1/s dt ≤ Cs Cs 2 s + 1 -1/s and ∞ Cs 2 s +1 Cs t 1/s dt = sCs 1 -s Cs 2 s + 1 1-1/s ≥ sCs 1 -s Cs 2 s + 1 -1/s , we have ∞ j= Cs 2 s +1 Csj -1/s ≤ 1 s ∞ Cs 2 s Cs t 1/s dt = C s s 1 -s 2-2s .Hence, we have1 n    Cs 2 s 2 + ∞ j= Cs 2 s +1 Csj -1/s    ≤ 1 n C s s + C s s 1 -s 2-2s . Set = 2 9/(4r+2s) κ 1/(2r+s) 2 -s 1 -s 1/(4r+2s) C s/(4r+2s) s n -1 4r+2s which is equivalent to 1 n C s s + C s s 1-s 2-2s = 1+2r 16 √2κ . Therefore, we attain the upper bound of ˜ n :˜ n ≤ 2 9/(4r+2s) κ 1/(2r+s) 2 -s 1 -s 1/(4r+2s) C s/(4r+2s) s n -1 4r+2s .Without loss of generality, we assume n is sufficiently large such thatn ≥ 2 9 κ 2 2 -s 1 -s C s s ⇔ 2 9/(4r+2s) κ 1/(2r+s) 2 -s 1 -s 1/(4r+2s) C s/(4r+2s) s n -1 4r+2s ≤ 1.Then n ≤ ˜ n ≤ 1. We now prove the following lemma. It is an extended version of Theorem 14.1 in [60] . Lemma A.6. We haveP    sup h H k ≤1 h 2 L 2 ρx,n -h 2 L 2 ρx h 2 L 2 ρx + t 2 ≤ 1 2    ≥ 1 -exp -c1nt 4r for any t ∈ [ n , 1]where c 1 is a constant independent of t n.Proof of Lemma A.6. We use a similar argument as in the proof of Theorem 14.1 in [60] . DefineZn(t) := sup h H k ≤1, h L 2 ρx ≤t h 2 L 2 ρx,n -h 2 L 2 ρx where t ∈ (0, 1]. Let E :=      sup h H k ≤1 h 2 L 2 ρx,n -h 2 L 2 ρx h 2 L 2 ρx + t 2 ≤ 1 2      c , A := Zn(t) ≥ t 2 2 , Ã := Zn(t) ≥ t 1+2r 2 .We first show that E ⊂ A. On the event E, there exists h ∈ H k such that h H k ≤ 1 andh 2 L 2 ρx,n -h 2 L 2 ρx h 2 L 2 ρx + t 2 > 1 2 . If h L 2 ρx ≤ t, we have h 2 L 2 ρx,n -h 2 L 2 ρx > 1 2 h 2 L 2 ρx + 1 2 t 2 ≥ 1 2 t 2 . Otherwise, set h = t h L 2 ρx h then h L 2 ρx = t and h 2 L 2 ρx,n -h 2 L 2 ρx > t 2 h 2 L 2 ρx • 1 2 h 2 L 2 ρx + 1 2 t 2 ≥ 1 2 t 2 . Therefore, E ⊂ A. Since t 2 2 ≥ t 1+2r 2for t ∈ (0, 1], we have E ⊂ A ⊂ Ã. To find an upper bound of EZ n (t), we use the symmetrization argument as follows:EZn(t) = E   sup h H k ≤1, h L 2 ρx ≤t 1 n n i=1 h(x i ) 2 -Eh(x) 2   = E   sup h H k ≤1, h L 2 ρx ≤t E 1 n n i=1 h(x i ) 2 - 1 n n i=1 h(x i ) 2 X   ≤ E   sup h H k ≤1, h L 2 ρx ≤t 1 n n i=1 wi h(x i ) 2 -h(x i ) 2   ≤ 2E   sup h H k ≤1, h L 2 ρx ≤t 1 n n i=1 wih(x i ) 2  where w 1 , • • • , w n are i.i.d. Rademacher variables. By Lemma A.15,2E   sup h H k ≤1, h L 2 ρx ≤t 1 n n i=1 wih(x i ) 2   ≤ 4κ • Q n (t).For t ≥ n , we have 4κEQUATION)Since1 n n i=1 sup h H k ≤1, h L 2 ρx ≤t E (h(x i ) 2 -Eh(x) 2 ) 2 ≤ sup h H k ≤1, h L 2 ρx ≤t E h(x) 4 ≤ κ 2 t 2 ,Lemma A. 16 givesP Zn(t) ≥ EZn(t) + u 1+2r 4 ≤ exp - 1 16 n 2 u 2+4r 4κ 2 • nEZn(t) + 2nκ 2 t 2 + 2 3 κ 2 • 1 4 nu 1+2r ≤ exp -c1n u 2+4r t 1+2r ∧ u 2+4r t 2 ∧ u 1+2rwhere c 1 = 1 96κ 2 . Putting u = t givesP( Ã) ≤ P Zn(t) ≥ EZn(t) + t 1+2r 4 ≤ exp -c1n t 1+2r ∧ t ≤ exp(-c1nt 4r ), i.e., P(E c ) ≥ P( Ãc ) ≥ 1 -exp(-c 1 nt 4r ).Let us return to our problem. Consider the eventE c =      sup h H k ≤1 h 2 L 2 ρx,n -h 2 L 2 ρx h 2 L 2 ρx + ˜ 2 n ≤ 1 2     .By Lemma A.6, we haveP(E c ) ≥ 1 -exp(-c 1 n˜ 4r n )where c 1 is a constant that does not depend on n. DefineT := min t ∈ N : 1 √ ηEt ≤ ˜ n .By the definition, we can easily obtain the upper bound of T asT < 1 + 1 ηE˜ 2 n . Since R(•) is non-decreasing, R 1 √ ηET ≤ R(˜ n) ≤ 1 c l Qn(˜ n)for some absolute constant c l where the second inequality follows from Lemma A.2. Note thath 2 L 2 ρx,n -h 2 L 2 ρx ≥ - 1 2 h 2 L 2 ρx + ˜ 2 n ⇒ h 2 L 2 ρx,n ≥ 1 2 h 2 L 2 ρx - 1 2 ˜ 2 nfor all h such that h H k ≤ 1 on the event E c . Thus,Qn(˜ n) = E   sup g H k ≤1, g L 2 ρx,n ≤˜ n 1 n n i=1 wig(x i ) X   ≤ E   sup g H k ≤1, g L 2 ρx ≤2˜ n 1 n n i=1 wig(x i ) X   on the event E c . Set F = g : X → R : g H k ≤ 1, g L 2 ρx ≤ 2˜ n .Then the ranges of functions in F are contained in [-κ, κ] and F = -F. By Lemma A.17 we haveE   sup g H k ≤1, g L 2 ρx ≤2˜ n 1 n n i=1 wig(x i ) X   ≤ 2E   sup g H k ≤1, g L 2 ρx ≤2˜ n 1 n n i=1 wig(x i )   + c1κ˜ 1+2r n with probability at least 1 -exp(-c 1 n˜ 1+2r n ) ≥ 1 -exp(-c 1 n˜ 4r n ). Hence, R 1 √ ηET ≤ 2 c l • Q n (2˜ n) + c1κ c l ˜ 1+2r n ≤ 1 κc l + c1κ c l ˜ 1+2rn holds with probability at least 1 -2 exp(-c 1 n˜ 4r n ). Here, the second inequality follows from (13) . Therefore,ER 1 √ ηET 2 ≤ 8 κc l + κc1 c l 2 ˜ 2+4r n • (0 ∨ (1 -2 exp(-c1n˜ 4r n ))) + κ 2 • 2 exp(-c1n˜ 4r n ) ≤ 8 κc l + κc1 c l 2 ˜ 2+4r n + 2κ 2 exp(-c1n˜ 4r n ) since R 1 √ ηET ≤ κ. From the fact that exp(-c1n˜ 4r n ) ˜ 2+4r n n 2r+1 2r+s exp(-c 1 n s 2r+s ) 1, we have ER 1 √ ηET 2 1/2 ˜ 1+2r n n -2r+1 4r+2s .Note that1 T ˜ 2 n n -1 2r+s and T ≤ 1 + 1 ηE˜ 2 n n 1 2r+s .Therefore, we bound the expected risk asE ιρ x (fT -f * 0 ) L 2 ρx B r-1/2 1 T + n -1 2r+s 1/2 n -r-1/2 2r+s + T 1/2 + n -1/2 2r+s T • ER 1 √ ηT E 2 1/2 + B r (1 + log T + T 1/2 n -1/2 2r+s )n -r 2r+s B r n -r 2r+s log n.A.3 Corollary of Theorem 3.4As mentioned in Section 3.3, one can remove B r in the upper bound in Theorem 3.4 by using more public inputs. The precise statement is as follows:Corollary A.7. Under Assumption 3.1, 3.2, and 3.3, with n 0 ≥ B 1+ n 1 2r+s (log(Bn)) 3 public inputs independently generated from ρx satisfying (2) DCL-KR gives the performance guaranteeE ι ρx (f j,T -f * 0 ) L 2 ρx ≤ C • n -r 2r+s log n for all j = 1, • • • , m where > 0 is a fixed constant, η ∈ (0, 1/κ 2) is a fixed learning rate, T is an adequate stopping rule, and the prefactor C does not depend on B, m, and n.Proof. In the proof of Theorem 3.4, there are two terms in the upper bound affected by B. One isB r-1/2 (log n 0 ) 3 n 0 r-1/2in the norm bound of the first term in (7) . The other isE (T k, ρx + λI) 1/2 (I -P Dp ) 2r ≤ 2 r λ r + (κ 2 + 1) • 4 exp - 1 4(κ 2 + 1)B 0in the norm bound of the third and fourth terms in (7) . For the first part,n 0 ≥ B 1+ n 1 2r+s (log(Bn)) 3 implies B r-1/2 (log n 0 ) 3 n 0 r-1/2 ≤ B -(r-1/2) n -r-1/2 2r+s (log(Bn)) -3r+3/2 (1 + ) log B + 1 2r + s log n + 3 log log(Bn) 3r-3/2 n -r-1/2 2r+s .For the latter part, set λ = 128(κ 2 + 1)2 n -1 2r+s /B. Then B 0 ≤ log κ 2 e + log(1/λ) λn 0 + log κ 2 e + log(1/λ) λn 0 ≤ 1 4(κ 2 + 1) log(Bn)and henceE (T k, ρx +λI) 1/2 (I -P Dp ) 2r ≤ 2 r λ r +(κ 2 +1)•4 exp - 1 4(κ 2 + 1)B 0 B -r n -r 2r+s + 1 Bn .Since it eliminates B r in the upper bound, we are done.Recall Cordes' inequality [16] . Lemma A.8 (Cordes' Inequality). Let A, B be two bounded positive linear operators on a seperable Hilbert space. Then for any s ∈ [0, 1] A s B s ≤ AB s holds.We also recall a property of projection operators. Lemma A.9 ( [53] ). Let Z be a bounded linear operator and P be a projection operator such that ran P = ran Z . Then for any bounded operator X and λ > 0 we have(I -P )X ≤ λ 1/2 (Z Z + λI) -1/2 X .There are some useful lemmas for PAC bounds. Lemma A. 10 ([17] ). Let X = {x 1 , • • • , x n } be a dataset where data points are independently generated from ν. Then(T k,ν + λI)(T k,X + λI) -1 ≤ 2 + 2 2κ 2 nλ + 4κ 2 N ν (λ) nλ log(2/δ) 2holds with confidence at least 1 -δ where δ ∈ (0, 1). Lemma A.11 ([38] ). Let X = {x 1 , • • • , x n } be a dataset where data points are independently generated from ν. Then(T k,ν + λI) -1 (T k,X + λI) ≤ 1 + 2κ 2 nλ + 4κ 2 N ν (λ) nλ log(2/δ)holds with confidence at least 1 -δ where δ ∈ (0, 1). Lemma A.12 ( [48] ). Let X = {x 1 , • • • , x n } be a dataset where data points are independently generated from ν. For λ ∈ (0, 1] such that N ν (λ) ≥ 1, (T k,ν + λI)(T k,X + λI) -1 ≤ 2 holds with confidence at least 1 -δ where4 exp   - 1 4(κ 2 + 1) • 1 + log N ν (λ) λn + 1 + log N ν (λ) λn -1   ≤ δ < 1.To prove Lemma A.2 in Appendix A.2.5, we introduce a concentration inequality for Lipschitz functions. Lemma A. 13 ([60] ). Let X 1 , • • • , X n be independent random variables whose supports are contained in [a, b] and f : R n → R be convex and L-Lipschitz with respect to the Euclidean norm. Then we haveP(f (X) ≥ Ef (X) + t) ≤ exp - t 2 4L 2 (b -a) 2 where X = [X 1 , • • • , X n ] and t > 0.Precisely, we use the following fact in Appendix A.2.5 Remark A.14 ([60] ). Let A ⊂ R n be a bounded set andf (x) = sup a∈A n k=1 a k x k where x = [x 1 , • • • , x n ] ∈ [-1, 1] n and a = [a 1 , • • • , a n ]. Since f (x) -f (x ) = sup a∈A n k=1 a k x k -sup a∈A n k=1 a k x k ≤ sup a∈A a, x -x R n ≤ sup a∈A a R n x -x R n ,f is a sup a∈A a R n -Lipshitz function where • R n is the Euclidean norm on R n . We can observe that f is convex since f is a supremum of convex functions defined on a convex compact set.To prove Lemma A.6 in Appendix A.2.5, we recall the Ledoux-Talagrand contraction inequality [29, 56] and Talagrand's inequality [4, 56] .Lemma A.15 (Ledoux-Talagrand Contraction Inequality). If φ : R → R is a L-Lipshitz function, then E sup h∈F 1 n n i=1 i φ(h(x i )) ≤ L • E sup h∈F 1 n n i=1 i h(x i ) .Lemma A.16 (Talagrand's Inequality). Let X 1 , • • • , X n be independent X -valued random variables. Let F be a countably family of measurable real-valued functions on X such that f ∞ ≤ U < ∞ and Ef (X i ) = 0 for all f ∈ F. LetZ := sup f ∈F n i=1 f (X i ), σ 2 ≥ 1 n n i=1 sup f ∈F E[f (X i ) 2 ], ν n := 2U EZ + nσ 2 .ThenP(Z ≥ EZ + t) ≤ exp - t 2 2ν n + 2 3 U t for all t ≥ 0.Lastly, we recall the following well-known property used in Appendix A.2.5. Lemma A.17 ([2] ). Let F be a class of functions with ranges in [a, b] and w 1 , • • • , w n be i.i.d. Rademacher variables. Then1 n E sup h∈F n i=1 w i f (x i ) ≤ inf α∈(0,1) 1 1 -α 1 n E sup h∈F n i=1 w i f (x i ) X + (b -a) log(1/δ) 4nα(1 -α)holds with probability at least 1 -δ. Also,1 n E sup h∈F n i=1 w i f (x i ) X ≤ inf α>0 (1 + α) 1 n E sup h∈F n i=1 w i f (x i ) + (b -a) log(1/δ) 2n 1 2α + 1 3holds with probability at least 1 -δ.As we mentioned before, DCL-NN considers the same problem as in Section 3 but local models are heterogeneous neural networks. That is, there are m parties andD i = {(x j i , y j i ) : j = 1, • • • , n i } is the private dataset of the ith party (i = 1, • • • , m) where D = m i=1 D i are i.i.d. whose distribution is ρ x,y .One remark is that the local data distributions of parties are not the same in general. To communicate training information, we introduce an unlabeled public input datasetZ = {z 1 , • • • , z n0 } ⊂ X .The goal of parties is to find a minimizer of the population risk E defined in Section 3.To extend DCL-KR to heterogeneous neural network settings, it is necessary to ensure that the assumptions of DCL-KR are satisfied as much as possible. Specifically, one important assumption in DCL-KR is the equality of kernels across local models. Indeed, public data predictions can vary in conflicting directions after the local training procedure, even when using the same local datasets, if the kernels differ. For further explanation of this claim, we consider the simple case of E = 1 in DCL-KR where E is the number of local iterations. After the consensus prediction u is distributed to local parties, the server then receives the updated local prediction on Z:EQUATIONfrom the ith local party. (The notation is consistent with Appendix A) Suppose two parties have exactly the same dataset. If the same kernel is used in these two parties, the updated local predictions f i (•) = w i g i (•) + b i (i = 1, • • • , m).3: # Feature Kernel Distillation Procedure 4: For party i (i = 1, • • • , m), compute the feature kernel values on public inputs {k fi (z j1 , z j2 ) :1 ≤ j 1 , j 2 ≤ n 0 } via (3) and upload them to the server. 5: The server aggregates the local feature kernel values to the target kernel values {k(z j1 , z j2 ) :1 ≤ j 1 , j 2 ≤ n 0 } with (4) and distributes them to all parties. 6:for party i = 1, • • • , m do 7: for t k = 0, • • • , T k -1 do 8: for mini-batch Z 0 ⊂ Z do 9:Update parameters of its model f i by maximizing CKA(k, k fi ) on Z 0 defined as (16) via gradient descent where k is fixed. for party i = 1, • • • , m do 17:Update w i and b i by maximizing MSE (Mean Squared Error) on the public inputs Z with consensus prediction y p,t via gradient descent with sufficiently many iterations. Upload the local prediction y i p,t+1 on Z to the server.end forThe server aggregates the local predictions to compute the consensus prediction y p,t+1 = m i=1 ni n y i p,t+1 and distributes y p,t+1 to all parties. 24: end for will be identical. However, if the kernels are different, this will not be the case. For kernels like the Gaussian kernel, which have high correlation between close inputs, the updated local predictions will be strongly influenced by data points close to each input. On the other hand, for kernels like the linear kernel, which have high correlation between distant inputs, the updated local prediction on Z will be influenced more by data points farther from each input. We can observe this fact from the above formula (14) . This observation implies that aggregating local learning information becomes very challenging when the kernels differ. In short, using the same kernel ensures that the shift mechanisms of predictions on Z at the edges are identical, making it possible for the aggregation through simple weighted averaging to work well. This is a key element of the strong theoretical results of DCL-KR and explains why kernel matching between neural networks is necessary in DCL-NN.Let f i be a local model of the ith party such thatf i (•) = w i g i (•) + b i , g i : X → R ci , w i ∈ R ci , c i ∈ N, and b i ∈ R for i = 1, • • • , m.Since most modern neural network architectures have a linear layer as the last layer, this setting is general enough. As (3), we set the feature kernel of f i(i = 1, • • • , m) to be k fi (x 1 , x 2 ) = g i (x 1 ) g i (x 2 ), x 1 , x 2 ∈ X .To bring the setting to the DCL-KR scheme, DCL-NN matches k f1 , • • • , k fm via kernel distillation procedure. Obviously, the target kernel in this procedure is a key factor in enhancing performance.Theoretically, using the ensemble kernelEQUATIONcan be a good way to construct a good kernel derived from local feature kernels. The reason is that this ensemble kernel is identical to the kernel induced by the (scaled) concatenation of the local feature maps, i.e.,k(x 1 , x 2 ) = m i=1 n i n k fi (x 1 , x 2 ) = m i=1 n i n g i (x 1 ) g i (x 2 ) = n1 n g 1 (x 1 ) n2 n g 2 (x 1 ) • • • nm n g m (x 1 )     n1 n g 1 (x 2 ) n2 n g 2 (x 2 ) • • • nm n g m (x 2 )     .In other words, the ensemble kernel has greater expressive power than individual feature kernels, and with a sufficient amount of data, it leads to better performance. We empirically verify that the performance of this ensemble kernel surpasses that of individual feature kernels in Figure 4 .DCL-NN sets this ensemble kernel k as the target kernel and local parties match their local feature kernels k f1 , • • • , k fm with the kernel k using the public dataset Z. For this purpose, we introduce Centered Kernel Alignment (CKA) [8] as a kernel similarity measure. The CKA between two kernels k 1 and k 2 on the public input distribution ρx is given byCKA(k 1 , k 2 ) = HSIC(k 1 , k 2 ) HSIC(k 1 , k 1 )HSIC(k 2 , k 2 )where HSIC(•, •) is a Hilbert-Schmidt Independence Criterion (HSIC) defined asHSIC(k i , k j ) = E x 1 ,x 2 ∼ ρx [k c i (x 1 , x 2 )k c j (x 1 , x 2 )] and the centered kernel k c i is given by k c i (x 1 , x 2 ) = k i (x 1 , x 2 ) -E x2 ∼ ρx [k i (x 1 , x2 )] -E x1 ∼ ρx [k i (x 1 , x 2 )] + E x1 ,x 2 ∼ ρx [k i (x 1 , x2 )], x 1 , x 2 ∈ X (i = 1, 2). However, since we have a finite number of samples, we employ the empirical CKA. The empirical CKA between two kernels k 1 and k 2 on inputs {c 1 , • • • , c p } is given byEQUATIONwhereK 1 =    k 1 (c 1 , c 1 ) • • • k 1 (c 1 , c p ) . . . . . . . . . k 1 (c p , c 1 ) • • • k 1 (c p , c p )    and K 2 =    k 2 (c 1 , c 1 ) • • • k 2 (c 1 , c p ) . . . . . . . . . k 2 (c p , c 1 ) • • • k 2 (c p , c p )   are Gram matrices and HSIC is an estimator of HSIC defined as While (empirical) CKA is a good metric for kernel matching, it is invariant to scaling, and therefore the local feature kernels resulting from the kernel distillation may have different scales. This affects the degree to which each local training influences during the DCL-KR-like follow-up procedure. To illustrate this point, consider the following example: Let the feature kernels of two local models f 1 and f 2 be as follows: k f1 (x, y) = φ(x) φ(y), k f2 (x, y) = (αφ(x)) (αφ(y)).HSIC(L, M ) = 1 (p -1) 2 tr(LHM H), L, M ∈ R p×pAfter distilling on the public data with consensus predictions, these two models become likef 1 (•) = w φ(•) + b, f 2 (•) = 1 α w (αφ(•)) + b = w φ(•) + b,i.e., two models are the same. Nevertheless, a gradient descent update on a data point (x 0 , y 0 ) with a learning rate η isw 1 ← w -η(w φ(x 0 ) + b -y 0 )for f 1 andw 2 ← 1 α w -η(w φ(x 0 ) + b -y 0 )for f 2 . Thus, after the gradient descent step we havef 1 (•) = (w -η(w φ(x 0 ) + b -y 0 )) φ(•) + b and f 2 (•) = 1 α w -η(w φ(x 0 ) + b -y 0 ) (αφ(•)) + b = (w -αη(w φ(x 0 ) + b -y 0 )) φ(•) + b.Hence the scale α affects the collaborative learning procedure. To address this issue, we compute the scale α using the estimator HSIC. Specifically, at the beginning of the collaborative learning phase, we compute α i = HSIC(K i , K i ) where K i is the Gram matrix with respect to the feature kernel of the ith party onZ (i = 1, • • • , m).Then we set the learning rate for the ith party asη 0 • max 1≤j≤m α 1/2 j α 1/2 iwhere η 0 is a base learning rate. In practice, computing HSIC over all public inputs is costly and unnecessary. Using only a small subset of public inputs is sufficient.We present DCL-NN in Algorithm 2. Here are some remarks.(1) The feature kernel distillation procedure requires only one round of two-way communication between the server and the parties.(2) The collaborative learning procedure follows the same process as in DCL-KR with g 1 , • • • , g m fixed. Note that kernel gradient descent reduces to standard gradient descent since the kernels have finite rank. In this process, if possible, optimization on the public dataset can be performed using the closed-form solution of kernel linear regression instead of gradient descent.(3) In this work, we apply FedMD for the pretraining of DCL-NN in the experiment. However, DCL-NN is a general algorithm that can use any algorithm for pretraining to obtain good local feature kernels. For example, kernel learning techniques [63] may be applied for pretraining.(4) Our algorithm can be naturally extended to regression problems with multi-dimensional outputs.For all datasets, we follow Algorithm 3 to construct non-i.i.d. settings. Note that this procedure is similar to the non-i.i.d. data generation procedure in classification tasks [44, 69] .Algorithm 3 Data Generating Procedure 1: Inputs: the number of parties m, the total number of private data n, the number of public inputs n 0 , the partition of input space A = {A 1 , • • • , A c } 2: Generate the whole private dataset size of n from ρ x,y and the public inputs size of n 0 from ρx . Sample two elements uniformly from the partition A. Set C ij = 1 if A i is chosen at the jth party and C ij = 0 otherwise. 9: end while 10:for i = 1, • • • , c do 11:Put the private data points where inputs are in A i into the datasets of parties with ratio Let X = [0, 1] ⊂ R and ρ x be the uniform distribution on X . The spaceα k C ik m j=1 α j C ij k=1,••• ,mH 1 := f ∈ AC[0, 1] f (0) = 0, f (x) 2 dρ x (x) < ∞is the reproducing kernel Hilbert space associated to the kernel k(x, y) = min(x, y) where AC[0, 1] is the collection of all absolutely continuous functions on [0, 1] [33, 60] . As mentioned in [33] , the covariance operator T k,ρx has eigenpairs {(λ i , e i )} i∈N whereλ i = 2i -1 2 π -2 , e i (x) = √ 2 sin 2i -1 2 πx .Thus, the eigenvalue decay rate s is 1 2 . Set a target functionf * 0 (x) = ∞ i=1 e i (x) i 3 = ∞ i=1 √ 2 i 3 sin 2i -1 2 πx .From the fact thatT 1/2-r k,ρx ∞ i=1 h i e i H k ≤ R ⇔ ∞ i=1 h 2 i λ 2r i ≤ R 2 , we have f * 0 = T 1/2 k,ρx g * 0 such that g * 0 H k = ∞ i=1 1 i 6 2i -1 2 π 4 1/2 =: R < ∞.Then r = 1. We generate data points from ρ x • ρ y|x such that ρ x is the uniform distribution on X as above and ρ y|x = N (y|f * 0 (x), 0.44 2 ).We divide X into a partition A := The experiments are implemented in PyTorch. We simulate a decentralized setting on a single deep learning workstation (Intel(R) Xeon(R) Gold 6430 with one NVIDIA GeForce RTX 4090 GPU and 189GB RAM). All DCL-KR implementations take less than 10 minutes per simulation. The non-parallel implementation of DCL-NN for large-scale datasets is completed within 48 hours. With parallel computing, the execution time for the same setup is expected to be reduced to within 2 hours.In this experiment, we evaluate DCL-KR by comparing its performance against two central models and three baselines. Specifically, we employ DC-NY [66] , DKRR-NY-CM [67] , and IED [48] as baselines. We also compare DCL-KR with central Kernel Ridge Regression (centralKRR) and central Kernel Regression with Gradient Descent (centralKRGD). As mentioned earlier, we evaluate the performance of these algorithms on Toy-1D and Toy-3D datasets.There are several hyperparameters for kernel machine-based algorithms: the ridge regularization hyperparameter λ in the ridge regressions and the number of iterations (or communication rounds) in the gradient descent-based regressions. We set λ= C • n -1 2r+s and T = int(D • n 1 2r+s) which are the optimal choices from theory. We determine the best values for C and D by grid search in our experiments (see Table 4 for the selected hyperparameter values). For centralKRGD and DCL-KR, we set the learning rate η = 0.5 which satisfies η ∈ (0, 1/κ 2 ) in Theorem 3.4. The number of local iterations E for DCL-KR is set to 5. For DKRR-NY-CM, we modify its Newton-Raphson iteration as shown below due to an instability issue, and we set the learning rate η = 0.01:u ← u -η • m j=1 n j n (P Z T k,Xi P Z + λI) -1 ((P Z T k,X P Z + λI)u -P Z S D y).The communication round T for DKRR-NY-CM is set to 10.To measure performance, we sample data from the distribution presented in Appendix C.1 for each simulation and compute E ι ρx (f i,T -f * 0 ) L 2 ρx by averaging the Root Mean Squared Errors (RMSEs) on the test dataset. We conduct 500 simulations for each setting, and the results are summarized in Figure 1 .As shown in Table 1 , DC-NY and DKRR-NY-CM have theoretical performance guarantees under the statistically homogeneous condition for a limited number of parties. However, they do not exhibit sufficiently good performance in massively distributed statistically heterogeneous settings. The performance degradation of DKRR-NY-CM appears to be linked to its second-order optimization scheme, which leads to ineffective batching in statistically heterogeneous settings. The performance degradation of DC-NY is expected, given the inherent limitations of divide-and-conquer algorithms.In contrast, IED demonstrates relatively better performance, despite the strong assumptions underlying its theory. Nevertheless, it still exhibits performance degradation compared to centralized models. DCL-KR, on the other hand, achieves performance comparable to centralized models, validating both the theoretical results and its practical feasibility.As public inputs directly affect the training information sharing, we anticipate that the performance of DCL-KR will vary depending on the number of public inputs n 0 . To examine this effect, we The results are summarized in Figure 2 . Consistent with theoretical results, the value of α does not affect the convergence rate of IED and DCL-KR. However, IED displays significant performance variations across different α values. In contrast, DCL-KR achieves its maximum performance when α is not too small (i.e., α ≥ 0.3 in Figure 2 ). This implies that DCL-KR requires fewer public inputs to achieve its maximal performance compared to IED, as predicted by theoretical results.So far, we consider the settings with ρ x = ρx . However, Theorem 3.4 covers the general case where ρx = ρ x . To verify this, we define the public input distribution in Toy-3D with the following density function (parametrized by β):p(x 1 , x 2 , x 3 |β) = 3 i=1 ((2 -2β)x i + β), (x 1 , x 2 , x 3 ) ∈ [0, 1] 3 , β ∈ (0, 1]. The Radon-Nikodym derivative dρx d ρx satisfies 0 ≤ dρx d ρx ≤ ( 1 β ) 3 .We conduct additional experiments to verify Theorem 3.4 and Corollary A.7, considering the case where β = 0.5 (with α = 1) and the case where β = 0.5 but α = 4 to compensate. The results are provided in Figure 5 . In the log-scale plot, the slope represents the convergence rate.First, we observe that the convergence rate of DCL-KR remains unchanged when β is changed from 1 (i.e., ρ x = ρx ) to 0.5. This observation is consistent with Theorem 3.4. Additionally, regarding Corollary A.7, we confirm that DCL-KR achieves performance almost identical to the case of ρ x = ρx by increasing n 0 . In contrast, the convergence rate of IED worsens when ρx changes, even when n 0 is increased. These experimental results highlight the advantages of DCL-KR in statistically heterogeneous environments. In the experiments on neural network-based collaborative learning algorithms, we evaluate three baselines (FedMD, FedHeNN, KT-pFL) and our algorithm DCL-NN. Note that while KT-pFL is a personalized collaborative learning algorithm, it also performs well in non-personalized settings, so we include it for comparison. Additionally, we evaluate centralized models as ideal cases and standalone models as worst cases. Centralized models are trained using all local data.We use two tabular datasets, Toy-3D and Energy, and three image datasets, RotatedMNIST, UTKFace, and IMDB-WIKI. The number of parties is set to 50 for all settings. For the tabular datasets, we employ four different fully connected neural networks (FNNs) with a ratio of 30%, 30%, 20%, and 20%. Specifically, There are fifteen 4-layer FNNs with 32 hidden units, fifteen 4-layer FNNs with 64 hidden units, ten 5-layer FNNs with 32 hidden units, and ten 3-layer FNNs with 64 hidden units. Similarly, for the image datasets, we use four different convolutional neural networks (CNNs) with the same ratio of 30%, 30%, 20%, and 20%. For the large-scale image datasets (RotatedMNIST and IMDB-WIKI), we construct 50 local parties using fifteen ResNet-18, fifteen ResNet-34, ten ResNet-50 [20] , and ten MobileNetv2 [54] . For UTKFace, which is an image dataset with limited data, we utilize four simpler CNN architectures due to the ineffectiveness of knowledge distillation with large underperforming models. The first and third CNNs share a similar architecture, featuring two convolutional layers with batch normalization [22] , two max pooling layers, and two fully-connected layers at the end. They differ only in the number of channels. In contrast, the second and fourth CNNs are more complex, featuring four convolutional layers with batch normalization, two max pooling layers, and two fully-connected layers at the end. They also differ in the number of channels. We use the ReLU activation function throughout. All optimizers used are Adam [26] . 2 One remark is that baseline algorithms only utilize a subset of public inputs through random sampling in each communication round, as performance tends to deteriorate due to overfitting when all public inputs are used in every round. Hyperparameters are tuned via grid search.Standalone models are trained with cross-validation and early stopping to prevent overfitting. To evaluate centralized models, we first compute the averaged test Root Mean Squared Error (RMSE) from at least 10 simulations for each neural network architecture and then calculate the weighted average of the performances of all architectures according to their ratio. For standalone models, we use the average of the test RMSEs of local models with the hyperparameters listed in Table 5 . In the table, RotatedMNIST is abbreviated as MNIST, and IMDB-WIKI is abbreviated as IMDB.For FedHeNN, we set the number of local epochs to 30 in all experiments. For KT-pFL, we set the number of local epochs to 10, the distillation coefficient to 0.5, and the learning rate of knowledge coefficient to 1e-3. Lastly, for DCL-NN in the main experiment, we set the learning rate for kernel matching to 1e-4 and the number of communication rounds in the collaborative learning phase to 50. We use the closed-form solution to train public data and full-batch gradient descent to train local data in the collaborative learning phase of DCL-NN and utilize FedMD for pretraining in DCL-NN.In the pretraining phase, the hyperparameters are the same as in the FedMD setting, except that the number of communication rounds is 100 for tabular data and 50 for image data. The remaining hyperparameters are presented in Table 6 , 7, 8, and 9. For all distillation-based collaborative learning algorithms, we simulate each setting at least 5 times with different initializations.Communication Efficiency Compared with FedMD and KT-pFL, DCL-NN incurs higher communication costs of O(n 2 0 ) due to the transmission of the Gram matrix. FedHeNN also utilizes kernel matching but performs it in batches for each communication round. Thus, in scenarios requiring many communication rounds, DCL-NN is more efficient than FedHeNN. However, pretraining also demands more communication cost. We leave the study of communication-efficient methods in DCL-NN for future work.In practice, the public inputs can be sampled from a distribution whose support is disjoint from that of the whole local input distribution. In this case, the assumption of DCL-KR does not hold; however, DCL-NN can still be applicable. To evaluate the performance of DCL-NN under these conditions, we compare the performance of DCL-NN and FedMD when the distribution of public inputs differs, as in Table 10 . We use CIFAR10 [28] for public inputs on UTKFace. As shown in Table 10 , DCL-NN Table 10 : Performance comparison of FedMD and DCL-NN on UTKFace with different public datasets. In addition to performance, the kernel performance of local feature kernels, computed in the same way as before, is shown in parantheses. does not yield better results in this case. Note that the kernel performance of local feature kernels is improved compared to FedMD, leading us to conclude that the performance degradation of DCL-NN with CIFAR10 is due to the violation of the DCL-KR assumption rather than the ineffectiveness of the kernel distillation procedure.Privacy Benefits of Distillation-based Collaborative Learning Due to its black-box nature, distillation-based information interaction is expected to offer privacy preservation benefits compared to parameter exchange (mainly done in FL) as mentioned in [19] . To the best of our knowledge, there is no rigorous study that discusses the privacy preservation advantages of distillation-based collaborative learning. We hope to see further discussion on this as well.Public Input Distribution Theorem 3.4 covers the case where the public input distribution ρx differs from that of local data inputs ρ x , but at least the support of ρx must include the support of ρ x . Therefore, we experimentally observe a performance drop in DCL-NN, a practical extension of DCL-KR, when ρx and ρ x have different supports. Enhancing the robustness of DCL-NN in this scenario is considered a promising direction for future work. Our theory does not cover situations where collecting public inputs is difficult. In such cases, a seperate generative model is usually trained to generate public inputs [72] . We leave the theoretical discussion that includes these cases for future work.• The authors should answer "Yes" if the results are accompanied by error bars, confidence intervals, or statistical significance tests, at least for the experiments that support the main claims of the paper. • The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall run with given experimental conditions). • The method for calculating the error bars should be explained (closed form formula, call to a library function, bootstrap, etc.) • The assumptions made should be given (e.g., Normally distributed errors).• It should be clear whether the error bar is the standard deviation or the standard error of the mean. • It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error bar than state that they have a 96% CI, if the hypothesis of Normality of errors is not verified. • For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are out of range (e.g. negative error rates). • If error bars are reported in tables or plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables in the text.• The answer NA means that the paper does not include experiments.• The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage. • The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute. • The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper).Justification: The authors verify that the research is conducted in the paper conform, in every respect, with the NeurIPS Code of Ethics.Guidelines:• The answer NA means that the authors have not reviewed the NeurIPS Code of Ethics.• If the authors answer No, they should explain the special circumstances that require a deviation from the Code of Ethics. • The authors should make sure to preserve anonymity (e.g., if there is a special consideration due to laws or regulations in their jurisdiction).Note that while DCL-NN should use vanilla gradient descent according to DCL-KR, Adam performs better in practice. 


Decision: The paper studies distillation-based collaborative learning (DCL) in a massively distributed heterogeneous setting. The authors analyze DCL kernel regression and prove its near-minimax optimality in their theoretical results. The paper also presents supporting experiments to demonstrate the effectiveness of the proposed method.

There is a consensus among reviewers that this work is theoretically sound and provides complementary experimental results to back the theoretical analysis. There are several suggestions by the reviewers to improve the paper further, especially the ones regarding the paper organization and the addition of discussions on (i) potential shortcomings if the local data belongs to a private domain without appropriate public data, and (ii) comparison of the theoretical bound in this work to those in prior works. 

Overall the reviews consistently agree that the paper makes a significant contribution over the prior work and provides a useful contribution to the domain. I would request the authors to consider the discussions summarized above and appropriately update the draft before the final submission.

Reviews: review: paper_summary: paper_summary: In this paper, the authors perform a study of distillation-based collaborative learning (DCL) in massively distributed statistically heterogeneous environments. In particular, the study focuses on analyzing DCL-KR, a non-parametric version of FedMD, and proves its near-minimax optimality. The authors proposed DCL-NN, a method that leverages kernel matching to implement DCL with heterogeneous neural networks. Experimental results support the theoretical analysis of DCL-KR and demonstrate the effectiveness of DCL-NN.
strength_weakness: strength_weakness: The paper does not perform any formal privacy analysis. As such, I believe using terms like "privacy-preserving" can be misleading in this context.
questions: questions: Could the authors avoid mentioning privacy as no formal privacy study is conducted?
limitations: limitations: The paper unnecessarily emphasizes privacy while admitting that privacy is not formally defined (cf. question). It would be nice to provide more discussion on this aspect or to present some empirical evidence of the claimed benefit (e.g., performing privacy attacks such as membership inference)
paper_summary: paper_summary: This paper theoretically proves a nonparametric version of the most standard distillation based collaborative learning algorithm (named DCL-KR) is nearly minimax optimal in massively distributed statistically heterogeneous environments.
Extensive experiments demonstrate their theoretical results and show the practical feasibility of DCL-NN.
strength_weakness: strength_weakness: Please see Limitations
questions: questions: Please see Limitations
limitations: limitations: This is not my field of study and I can only give a Borderline Reject. I will adjust my score based on the scores of other more specialised judges.
paper_summary: paper_summary: The paper investigates the theoretical underpinnings and practical implementation of distillation-based collaborative learning (DCL) from a kernel regression perspective. The authors propose DCL-KR, a nonparametric version of the FedMD algorithm, which achieves nearly minimax optimal convergence rates in massively distributed and statistically heterogeneous environments. Building on these theoretical insights, the authors introduce DCL-NN, a practical DCL algorithm designed for neural networks, which incorporates feature kernel matching to align local models. Extensive experiments on various regression tasks demonstrate the superiority of DCL-KR and DCL-NN over existing methods.
strength_weakness: strength_weakness: 1. The paper lacks a discussion on the privacy risk and communication efficiency of DCL-KR and DCL-NN.

2. The paper is based on distillation-based collaborative learning method that requires public data, which limits the applications on private domains.
questions: questions: 1. How does the choice of kernel affect the performance/theoretical guarantee of the proposed approach?

2. Is transferring predictions a common practice in collaborative learning for regression? If not, how does it compare with the other approaches in terms of privacy and efficiency?
limitations: limitations: Yes.
paper_summary: paper_summary: The authors propose a nonparametric version of FedMD, a distillation-based collaborative learning methodology. They also propose a neural network variant of the nonparametric approach as an extension for heterogeneous local neural network models. They provide both theoretical results on the nonparametric approach and experimental results demonstrating efficacy of the neural network extension. The theoretical results prove the nearly minimax optimality of this nonparametric collaborative learning algorithm in massively distributed statistically heterogeneous environments.
strength_weakness: strength_weakness: (1) The need to perform learning rate scaling to ensure the impact of local iterations is consistent across models introduces additional complexity. The exact method for scaling may vary and might not be trivial to determine or implement. (2) The computation of Gram matrices for CKA involves O(p^2) operations, where p is the number of public data points. This quadratic complexity can become prohibitive as the size of the public dataset increases. (3) The method assumes that the ensemble kernel k = \sum_{i=1}^{m} \frac{n_i}{n} k_{f_i} will perform better than individual feature kernels. This assumption might not hold in cases where local models are highly heterogeneous, as the averaged kernel could lose important characteristics of individual models, leading to suboptimal performance. The authors empirically verify the superiority of the ensemble in the appendix.
questions: questions: (1) How do you handle the computational burden associated with the calculation of Gram matrices for large datasets? (2) How do you handle cases where these representations before the last layers of the neural networks are suitable for kernel matching?
limitations: limitations: N/A
paper_summary: paper_summary: This work presents an analysis of a distillation-based collaborative learning algorithm called FedMD from a nonparametric perspective. The method is to adopt an operator theoretic approach to obtain an upper rate for the expected generalization error of local models. Then, the authors propose DCL-KR to achieve a privacy-preserving nearly minimax optimal collaborative learning algorithm with kernel regression in massively distributed statistically heterogeneous environments.
strength_weakness: strength_weakness: - The paper's organization and structure are not clear.  
- The writing looks like combining different parts of different papers without a clear clue.
- The paper suddenly gives many assumptions in Section 3 for analysis but does not give any intuition or insight into them.
- Some important experimental results are in the appendix. I would suggest that the authors extract insight into the experimental findings in the main body.
- The language used in the paper is not concise.
questions: questions: The related work only gives a rough summary of the previous works. When the authors focus on the theoretical analysis, could the authors compare their theoretical bounds to other works?

There are many important experiments in the appendix, such as C.3.1. Can the authors summarize empirical findings with a better paper structure? 

In Section 3.3, the paper starts with "To derive theoretical results" and lists a lot of assumptions without any intuitive understanding behind them. Could the authors give an overview and intuition of theoretical results?

Minor: I feel very confused when reading the paper as it introduces many conclusions/results without too much explanation. I suggest the authors elaborate on the intuitions before going into the details.
limitations: limitations: The paper presents a discussion of limitations, which is not convincing. For example, "there is no rigorous study that discusses the privacy preservation advantages of distillation-based collaborative learning". Actually, the standard technique for privacy is not distillation.

Comments: title: Thank you for the responses
comment: Thanks for the response. This has cleared up most concerns for me.

On this comment, 
>Therefore, with enough data, it is not possible for the ensemble kernel to perform worse than individual local feature kernels. However, fully distilling the ensemble kernel is challenging. Thus, there is a possibility of losing important characteristics (that the reviewer concerned with) during the kernel distillation process. Precisely, some important features may have a lower portion in the ensemble kernel and so they can be potentially missed during the distillation process. One possible solution is for each local party to construct local feature kernels by selecting important features and adjusting weights accordingly.

This is particularly what I was aiming at. On re-reading my review, I realise I was not clear.
I agree with enough data the ensembled kernel cannot be worse than the local kernels. I agree the problem here is the distillation process and that piece should be explored in future work. There is work in leveraging bochner's theorem for learning kernels from data in the gaussian process literature. Perhaps these approaches could provide more optimal kernels at the local level.

I've updated my score given this response.
After reading feedback from other reviewers, I am willing to improve my rating based on trusting the authors to modify the paper.
Thank you for the responses. I will keep my positive score.
Thanks for the detailed explanation and response. After reading the whole page, I am willing to improve my rating based on trusting the authors to modify the paper.

Thank you very much.
Thank you for the author's response. I will take into account the feedback from other reviewers to adjust my scoring subsequently.
Thanks for the author's response. The response has addressed my questions and I'll keep my positive score.