Year: 2021

Abstract:
In this paper, we study the statistical properties of distributed kernel ridge regression together with random features (DKRR-RF), and obtain optimal generalization bounds under the basic setting, which can substantially relax the restriction on the number of local machines in the existing state-of-art bounds. Specifically, we first show that the simple combination of divide-and-conquer technique and random features can achieve the same statistical accuracy as the exact KRR in expectation requiring only O(|D|) memory and O(|D|) time. Then, beyond the generalization bounds in expectation that demonstrate the average information for multiple trails, we derive generalization bounds in probability to capture the learning performance for a single trail. Finally, we propose an effective communication strategy to further improve the performance of DKRR-RF, and validate the theoretical bounds via numerical experiments.

heading: 1 INTRODUCTION
text: Kernel ridge regression (KRR) is one of the most popular nonparametric learning methods (Vapnik, 2000). Despite the excellent theoretical guarantees, KRR does not scale well in large scale settings because of high time and memory complexities (Liu et al., 2013; 2014; 2017; 2018; 2020b; Liu & Liao, 2015; Li et al., 2018; 2019c). Distributed learning (Zhang et al., 2013; Hsieh et al., 2014; Chang et al., 2017b; Li et al., 2019b; Lin et al., 2020), random features (Rahimi & Recht, 2007; Sutherland & Schneider, 2015; Rudi & Rosasco, 2017; Rudi et al., 2018; Liu et al., 2020a; Avron et al., 2017a; Yu et al., 2016; Jacot et al., 2020), and Nyström methods (Drineas & Mahoney, 2005; Ding & Liao, 2012; Yang et al., 2012; Camoriano et al., 2016; Si et al., 2016; Musco & Musco, 2017; Kriukova et al., 2017) are the most widely used large scale techniques to address the scalability issues. Recent statistical learning works on KRR together with large scale approaches demonstrate that these large scale approaches can not only obtain great computational gains but also can guarantee the optimal theoretical properties, such as KRR with divide-and-conquer (Zhang et al., 2013; 2015; Chang et al., 2017b;a; Guo et al., 2017; Lin et al., 2017; Li et al., 2019b;d; Lin et al., 2020), with random features (Rudi & Rosasco, 2017; Li et al., 2019e; Carratino et al., 2018; Yang et al., 2012), and with Nyström methods (Bach, 2013; Alaoui & Mahoney, 2015; Rudi et al., 2015; 2017; Ding et al., 2020).
The combinations of distributed learning and other large scale approaches are very intuitive but effective strategies to further improve the effectiveness, such as distributed learning with gradient descent algorithms (Lin & Zhou, 2018; Richards et al., 2020), with multi-pass SGD (Lin & Rosasco, 2017; Lin & Cevher, 2018; 2020), with random features (Li et al., 2019b), and with Nyström methods (Yin et al., 2020). The optimal generalization performance of these combining approaches has been studied, however, the main theoretical problem is that there is a strict restriction on the number of local machines. For sample, in (Lin & Zhou, 2018; Li et al., 2019b; Yin et al., 2020), to guarantee the optimal generalization performance in the basic setting, the upper bounds of the local machines are restricted to be a constant, which is difficult to be satisfied in real applications. ∗Both authors contributed equally to this work †Corresponding author
In this paper, we aim at enlarging the number of local machines by considering communications among different local machines. This paper makes the following three main contributions. Firstly, we improve the existing state-of-art results of the divide-and-conquer technique together with random features. We prove that the optimal generalization performance can be guaranteed even the partitions reach Ω( √ |D|), which are limited to a constant Ω(1) for the existing bounds in the basic setting, |D| is the size of the data sets. Secondly, to essentially reflect the generalization performance, beyond the minimax optimal rates in expectation, we derive optimal learning rates in probability, which can capture the learning performance for a single trail. Finally, we develop a communication strategy to further improve the performance of our proposed method, and validate the effectiveness of the proposed communications via both theoretical assessments and numerical experiments.
Related Work The most related work includes the statistical analysis of distributed learning and random features.
Distributed learning. Optimal learning rates for divide-and-conquer KRR in expectation were established in the seminal work (Zhang et al., 2013; 2015). An improved bound was derived in (Lin et al., 2017) based on a novel tool of integral operator. Based on the proof techniques proposed in (Zhang et al., 2013; 2015; Lin et al., 2017), optimal learning rates were established for distributed spectral algorithms (Guo et al., 2017; Blanchard & Mücke, 2018; Lin & Cevher, 2020), distributed gradient descent algorithms (Lin & Zhou, 2018; Richards et al., 2020), distributed semi-supervised learning (Chang et al., 2017b), distributed local average regression (Chang et al., 2017a; Lin & Cevher, 2020), localized SVM (Meister & Steinwart, 2016), etc. Some other communication strategies for distributed learning have been provided, see e.g. (Fan et al., 2019; Li et al., 2019a; Lin & Cevher, 2020; Li et al., 2020), and references therein. The theoretical analysis mentioned above shows that the divide-and-conquer learning can achieve the same statistical accuracy as the exact KRR, however, there is a strict restriction on the number of local machines. The optimal learning rates with a less strict condition on the number of local machines for distributed stochastic gradient methods and spectral algorithms were established in (Lin & Cevher, 2020). In (Lin et al., 2020), they considered the communications among different local machines to enlarge the number of local machines. However, the communication strategy proposed in (Lin et al., 2020) based on an operator representation, which requires communicating the input data between each local machine. Thus, it is difficult to protect the data privacy of each local machine. Furthermore, for each iteration, the communication complexity of each local machine is O(|D|d), where d denotes the dimension, which is infeasible in practice for large scale data sets.
Random Features. The generalization bound of random features was first proposed in (Rahimi & Recht, 2008), which shows that O(|D|) random features are needed for O(1/ √ |D|) learning rate. Some works further studied its theoretical performance (Cortes et al., 2010; Yang et al., 2012; Sutherland & Schneider, 2015; Sriperumbudur & Szabó, 2015; Bach, 2017; Avron et al., 2017b). By applying the standard integral operator framework (Smale & Zhou, 2007; Caponnetto & Vito, 2007), the optimal generalization bounds of KRR with random features were established in (Rudi & Rosasco, 2017), which requires only O( √ |D| log(|D|)) random features. To decrease the size of random features, an improved approach was proposed based on a novel leverage score sampling strategy (Rudi et al., 2018). Sun et al. (2018) extended the result of (Rudi & Rosasco, 2017) to SVM. In (Li et al., 2019e), they further devised a simple framework for the unified analysis of random Fourier features, which can be applied to KRR, as well as SVM and logistic regression. To further improve the effectiveness, recently, Li et al. (2019b) considered the simple combination of divide-and-conquer and random features. However, to guarantee the optimal generalization performance, the number of local machines should be restricted to a constant, degenerating it into a single random features-based large scale KRR.

heading: 2 BACKGROUND
text: In a standard framework of supervised learning, there is a probability space X × Y with an unknown distribution ρ, where X is the input and Y is the output space. The sample set D = {(xi, yi)}ni=1 of size n is drawn i.i.d from X × Y with respect to ρ. Let K : X × X → R be a Mercer kernel, and HK be its reproducing kernel Hilbert space (RKHS) (Steinwart & Christmann, 2008; Vapnik, 2000), and assume that K(x,x′) ≤ κ,∀x,x′ ∈ X . Throughout, we will denote the inner product inHK by 〈·, ·〉K , and corresponding norm by ‖ · ‖K .

heading: Kernel Ridge Regression (KRR)
text: KRR is one of the most popular nonparametric learning methods (Shawe-Taylor & Cristianini, 2000; Vapnik, 2000), which can be stated as
fD,λ = arg min f∈HK  1|D| |D|∑ i=1 (f(xi)− yi)2 + λ‖f‖2K  , (1) where λ > 0 is the regularization parameter, |D| is the size of D. Using the representation theorem (Shawe-Taylor & Cristianini, 2000; Vapnik, 2000), fD,λ can be written as fD,λ(x) =∑|D| i=1 αiK(xi,x) with α = (KD + λI) −1yD, where KD = 1|D| [KD(xi,xj)] |D| i,j=1 is the |D|× |D| kernel matrix and yD = 1|D| (y1, . . . , y|D|) T.
Despite the excellent theoretical guarantees (Blanchard & Krämer, 2010; Caponnetto & Vito, 2007), KRR requires O(|D|2) memory to store KD, and O(|D|3) time to solve inverse of KD + λI, which is infeasible for large scale settings.

heading: KRR with Random Features (KRR-RF)
text: Assuming the spectral measure has a density function %(·), the corresponding shift-invariant kernel can be written as K(x,x′) = ∫ Ω ψ(x,ω)ψ(x′,ω)%(ω)dω, where ψ : X × Ω → R is a continuous and bounded function with respect to ω and x. The main idea behind random Fourier features is to approximate the kernel function K(x,x′) by its Monte-Carlo estimation (Rahimi & Recht, 2007): KM (x,x′) = 1M ∑M i=1 ψ(x,ωi)ψ(x
′,ωi) = 〈φM (x),φM (x′)〉, where φM (x) =
1√ M (ψ(x,ω1), . . . , ψ(x,ωM )) T. The solution of KRR with random features can be
written as
fM,D,λ(x) = w T M,D,λφM (x) with wM,D,λ = (ΦM,DΦ T M,D + λI) −1ΦM,DȳD, (2)
where ΦM,D = 1√|D| (φM (x1), . . . ,φM (x|D|)) and ȳD = 1√ |D| (y1, . . . , y|D|) T.
KRR-RF requires O(M |D|) to store ΦM,D, O(M3) and O(M2|D|) time to solve the inverse of (ΦM,DΦ T M,D + λI) and the matrix multiplication ΦM,DΦM,D, respectively. Thus, the total space and time complexity of KRR-RF are O(M |D|) and O(M2|D|), M |D|, respectively.

heading: Distributed KRR with Random Features (DKRR-RF)
text: Let {Dj}mj=1 be m disjoint subsets with D = ∪mj=1Dj . The distributed KRR with random features (DKRR-RF) is defined as
f̄0M,D,λ = m∑ j=1 |Dj | |D| fM,Dj ,λ, (3)
where fM,Dj ,λ(x) = w T M,Dj ,λφM (x) with wM,Dj ,λ = (ΦM,DjΦ T M,Dj + λI) −1ΦM,Dj ȳDj . The space complexity, time complexity and communication complexity of DKRR-RF for each local machine are O(M |Dj |), O(M3 +M2|Dj |) and O(M), respectively.

heading: 3 DKRR-RF WITH COMMUNICATIONS (DRKK-RF-CM)
text: In this section, we will present an effective communication strategy to enlarge the number of local machines. We first give the motivation of our communication strategy, and then propose a communication-based method, called DRKK-RF-CM. The proposed communication strategy are adaptations from (Lin et al., 2020) to avoid communicating local data among partition nodes.
Motivation. Let gM,D,λ : RM → RM be gM,D,λ(w) := [ ΦM,DΦ T M,D + λI ] w −ΦM,DȳD. One
can see that 2gM,D,λ(w) is the gradient of the empirical risk of 1|D| ∑ (xi,yi)∈D ( wTφM (xi)− yi )2 + λ‖w‖2 on w. From Eq.2, we know that for any w, the following equation holds:
wM,D,λ = w − [ ΦM,DΦ T M,D + λI ]−1 [ [ ΦM,DΦ T M,D + λI ] w −ΦM,DȳD ] = w − [ ΦM,DΦ T M,D + λI ]−1 gM,D,λ(w).
(4)
Algorithm 1 Distributed KRR with Random Features and Communications (DKRR-RF-CM) Initialize: w̄0M,D,λ = 0 for t = 1 to p do
Local machine: compute the local gradient gM,Dj ,λ(w̄ t−1 M,D,λ), and communicate back to GM. Global machine: get the global gradient gM,D,λ(w̄t−1M,D,λ) = ∑m j=1 |Dj | |D| gM,Dj ,λ(w̄ t−1 M,D,λ) and communicate to each local machine. Local machine: compute βt−1j = [ ΦM,DjΦ T M,Dj + λI ]−1 gM,D,λ(w̄ t−1 M,D,λ) and communicate back to the global machine. Global machine: compute w̄tM,D,λ = w̄ t−1 M,D,λ − |Dj | |D| ∑m j=1 β t−1 j , and communicate to each
local machine. end for Output: w̄pM,D,λ and f̄ p M,D,λ = 〈w̄ p M,D,λ,φM (·)〉
Define w̄0M,D,λ = ∑m j=1 |Dj | |D| wM,Dj ,λ, it is easy to verify that
w̄0M,D,λ = w − m∑ j=1 |Dj | |D| [ ΦM,DjΦ T M,Dj + λI ]−1 gM,Dj ,λ(w). (5)
Comparing 4 and 5, and noting that the global gradient gM,D,λ(w) can be achieved via the communications of each local gradient gM,Dj ,λ(w), i.e., gM,D,λ(w) = ∑m j=1 |Dj | |D| gM,Dj ,λ(w), thus, we consider the following Newton Raphson iteration-based communication strategy:
w̄tM,D,λ = w̄ t−1 M,D,λ − m∑ j=1 |Dj | |D| [ ΦM,DjΦ T M,Dj + λI ]−1 gM,D,λ(w̄ t−1 M,D,λ). (6)
We propose an iterative procedure to implement the communication strategy Eq.6, which can be broken down into 4 steps. At first, each local machine computes the local gradient and communicates back to the global machine. Then the global machine computes the global gradient based on the local gradient, and communicates it to each local machine. In the third step, each local machine computes βt−1j , and communicates back to the global machine. Finally, the global machine obtains the solution w̄tM,D,λ. More details can be seen in Algorithm 1.
Complexity analysis. Space complexity: each local machine only needs to store ΦM,Dj and the local gradient gM,Dj ,λ, thus the space complexity of each local machine isO(M |Dj |+M) = O(M |Dj |); Time complexity: for each local machine, we only need to compute the matrix multiplication ΦM,DjΦ T M,Dj and the inverse of ΦM,DjΦ T M,Dj + λI once. For each iteration, we need to compute the local gradient gM,Dj ,λ and βj for each local machine. Therefore, the total time complexity of each local machine is O(M3 + M2|Dj | + pM |Dj |), where p is the number of communication; Communication complexity: for each iteration, we only communicate the local gradient gM,Dj ,λ and βj to the global machine, and receive the gradient gM,D,λ and w̄t−1M,D,λ from the global machine, so the total communication complexity is O(pM). Remark 1. From the complexity analysis above, we can see that if the number of the communication p satisfying p ≤M or p ≤ |Dj |, then the time and space complexity of DKRR-RF-CM are the same as DKRR-RF. Only the communication complexity is slightly increased from O(M) to O(pM).

heading: 4 THEORETICAL ANALYSIS
text: In this section, we analyze the generalization performances of DKRR-RF and DKRR-RF-CM. The performance of the algorithm is usually measured by the expected risk E(f) = ∫ X×Y(f(x) − y)2dρ(x, y). The optimal hypothesis fHK inHK is denoted by fHK = arg minf∈HK E(f), and we assume fHK exists in the paper.

heading: 4.1 OPTIMAL LEARNING RATE FOR DKRR IN EXPECTATION
text: Theorem 1. Suppose that |ψ(x,ω)| ≤ τ almost surely, τ ∈ [1,∞) and |y| ≤ ζ. If λ = Ω(|D|− 12 ), |D1| = . . . |Dm|, the number of partitions m and the number of random features M respectively correspond to m . √ |D| and M & √ |D|, then, for every δ ∈ (0, 1], with probability at least 1− δ,
we have E[E(f̄0M,D,λ)]− E(fHK ) = O ( |D|− 12 log2(1/δ) ) . From Theorem 1, one can see that if m . √ |D| and M & √ |D|, the learning rate of the general-
ization bound can reach O(1/ √ |D|), which is optimal in a minmax sense (Rudi & Rosasco, 2017; Caponnetto & Vito, 2007). It means that, in this basic setting, as long as the number of partitions and random features are in order Ω( √ |D|), the corresponding ridge regression estimator has optimal generalization properties. The assumption of |y| ≤ ζ can be related to the Bernstein condition or moment assumption (Blanchard & Krämer, 2016), but for simplicity, in this paper, we only consider that Y is bounded. Optimal learning rates for divide-and-conquer KRR in expectation have been established in (Zhang et al., 2013; 2015; Lin et al., 2017), etc. However, there is a strict restriction on the number of local machines m. Specifically, in (Lin et al., 2017), to reach the optimal rate, m should to restrict to a constant m = Ω(1). In (Li et al., 2019b), the authors have studied the generalization performance of the combination of divide-and-conquer technique and random features. Using the same setting as Theorem 1 (that is r = 1/2 and γ = 1 in Theorem 8 of Li et al. (2019b)), they prove that, if M & √ |D| and m . Ω(1), then E[E(f̄0M,D,λ)]− E(fHK ) = O(|D|− 1 2 log2(1/δ)). It means that to guarantee the optimal generalization properties, the number of partitions should be restricted to a constant, but for our result is Ω( √ |D|). In (Li et al., 2019b), they also considered using the unlabeled data to enlarge the number of partitions. They have proved that (see Corollary 12 of Li et al. (2019b) for detail), if M & √ |D| and m . |D∗|/|D|, then E[E(f̄0M,D,λ)]−E(fHK ) = O(|D|− 1 2 log2(1/δ)), where D∗ is the dataset includes both labeled and unlabeled data. Thus, if we want the number of partitions m to reach Ω( √ |D|) as the same as our Theorem 1, the size of |D∗| should be Ω( √ |D|)
times of |D|. In this case, the data size of each local machine is |D| = |D|3/2/ √ |D|, so the time and space complexity are the same as KRR with a single random feature technique.

heading: 4.2 OPTIMAL LEARNING RATES FOR DKRR IN PROBABILITY
text: Note that E[E(f̄0M,D,λ)]−E(fHK ) = E[‖f̄0M,D,λ−fHK‖2ρ] (Caponnetto & Vito, 2007), thus Theorem 1 proposes the optimal learning rate in expectation, which demonstrates the average information for multiple trails, but may fail to capture the learning performance for a single trail. To essentially reflect the generalization performance for a single trail, we derive the optimal learning rate in probability:
Theorem 2. Under the same assumptions as Theorem 1. If λ = Ω(|D|− 12 ), |D1| = . . . |Dm|, m . |D| 14 and M & |D| 12 , then, for every δ ∈ (0, 1], with probability at least 1 − δ, we have∥∥∥f̄0M,D,λ − fHK∥∥∥2
ρ = O
( |D|− 12 log2(1/δ) ) .
To guarantee the optimal generalization properties in probability, the number of partitions should be restricted to Ω(|D|1/4), which is stricter than Ω(|D|1/2) in Theorem 1. This is because the generalization error in expectation can be decomposed into approximation error, sample error and distributed error (more details can be seen in Proposition 1 in Appendix), but the error decomposition in probability is not easy to separate a distributed error in probability to control the number of local machines. The derive the optimal learning rate, we provide a novel decomposition, please see details in Proposition 9 in Appendix.
The following result demonstrates that the proposed communication strategy can enlarge the number of partitions in probability.
Theorem 3. Under the same assumptions as Theorem 1, If λ = Ω(|D|− 12 ), |D1| = . . . |Dm|, m . |D| p+1 2(p+2) and M & |D| 12 , then, for every δ ∈ (0, 1], with probability at least 1 − δ, we have
Table 1: Computational complexity required by different algorithms for the optimal learning rate
O
(
1/
√
|D|)
)
in the basic setting. Logarithmic terms are not showed.
Methods Partitions m Random M Types Space Time Communication
KRR (Caponnetto & Vito, 2007) / / In probability |D|2 |D|3 / KRR-RF (Rudi & Rosasco, 2017) / |D|0.5 In probability |D|1.5 |D|2 / KRR-Nystöm Rudi et al. (2015) / |D|0.5 In probability |D|1.5 |D|2 / DKRR (Zhang et al., 2015; Chang et al., 2017b) |D|0.5 / In expectation |D| |D|2 |D|0.5 DKRR (Lin et al., 2020) |D|0.25 / In probability |D|1.5 |D|2.25 |D|0.75 DKRR-CM (Lin et al., 2020) |D| p+1 2(p+2) / In probability |D| p+3 p+2 |D| 3(p+3) 2(p+2) pd|D| DKRR-RF (Li et al., 2019b) Ω(1) |D|0.5 In expectation |D| |D|2 |D|0.5 DKRR-RF (Theorem 1) |D|0.5 |D|0.5 In expectation |D| |D|1.5 |D|0.5 DKRR-RF (Theorem 2) |D|0.25 |D|0.5 In probability |D|1.25 |D|1.75 |D|0.5 DKRR-RF-CM (Theorem 3) |D| p+1 2(p+2) |D|0.5 In probability |D| 2p+5 2p+4 |D| 3p+7 2p+4 p|D|0.5
∥∥∥f̄pM,D,λ − fHK∥∥∥2 ρ = O ( |D|− 12 logp+2(1/δ) ) , where f̄pM,D,λ is returned by Algorithm 1 under p-th iterations.
Compared Theorem 3 with Theorem 2, it is clear that the proposed communication strategy can relax the restriction on m from Ω ( |D|1/4 ) to Ω ( |D|(p+1)/(2(p+2)) ) . Note that m is monotonically increasing with the number of communications p, which can demonstrate the power of the proposed communications. When p → ∞, the partitions can reach Ω( √ |D|), which is the same as the generalization bound in expectation.
Remark 2. In the main text of this paper, we only give the optimal rates of DKRR-RF in the basic setting. The fast learning rates can be achieved under favorable conditions, see in Appendix.
Remark 3 (The Significance of Distributed Learning for RF). At first glance, it seems that the bottleneck in learning with random Fourier features is not the size of the dataset but the number of features. However, from (Rudi & Rosasco, 2017), one can see that we only requiring O( √ |D|) random features to guarantee the optimal performance, so the total computational complexity is O(M3 + M2|D|) = O(|D|2). Thus, the computational bottleneck in learning is not only the size of random features, but also the size of dataset. If we don’t consider reducing the size of D, the computational complexity is |D|2 in the basic setting, which is not suitable for large scale problems. Distributed learning is one of the most popular methods to reduce the size of dataset. The distributed learning bring the distributed error, but can decrease the variance of the model Zhang et al. (2013; 2015). Thus, how to choose an appropriate number of partitions to trade off the distributed error and the variance to guarantee the optimal performance is a very interesting and significant direction.

heading: 4.3 COMPARED WITH THE RELATED WORK

heading: Comparisons of the Time and Space Complexities
text: Table 1 reports the statistical and computational properties of the related approaches and our theoretical findings under the basic setting. We see that our DKRR-RF can guarantee the optimal generalization performance in expectation only requiring |D| memory and |D|1.5 time, which is more effective than other methods. For DKRR-RF-CM, we can also see that it can guarantee optimal generalization performance in probability requiring less complexity than the communication-based method of DKRR-CM (Lin et al., 2020).
Remark 4. In (Rudi et al., 2017), the authors considered combining the Nyström method and preconditioned conjugate gradient (PCG) (Cutajar et al., 2016) to scale up KRR. As far as we know, it is the only existing work that can guarantee optimal statistical accuracy, only requiring |D| memory and |D|1.5 time for KRR. In this paper, we consider combining distributed learning and random features, a completely different path from (Rudi et al., 2017). Note that in our proposed method, we need to compute the inverse of ΦM,DjΦ T M,Dj + λI, which requires |D|
1.5 time. Inspired by (Rudi et al., 2017), we can also adopt PCG to avoid the inverse calculation, which can further speed up our proposed DKRR-RF-CM. The combination of DKRR-RF-CM and PCG may open a path to reach the linear time complexity for optimal learning rate.

heading: Novelty and Proof Techniques
text: The most related works of our paper are (Li et al., 2019b), (Lin et al., 2020) and (Rudi & Rosasco, 2017). We discuss the novel techniques adopted to derive the improved results compared with them.
Compared with (Li et al., 2019b). (a) To derive the learning bounds, QM,D := ‖(CM + λI)−1/2(CM − CM,D)(CM + λI)−1/2‖ is required to be estimated, where CM and CM,D are self-adjoint and positive operators defined in Definition 1 (see in Appendix). In (Li et al., 2019b), they used a classical approach from (Chang et al., 2017b; Guo et al., 2017) to estimateQM,D (see Lemmas 21 and 22 in (Li et al., 2019b)), and obtain that QM,D ≤ 1√λ‖(CM − CM,D)(CM + λI)
−1/2‖ = O ( 1/λ|D|+ √ N (λ)/λ|D| ) , where N (λ) is the effective dimension defined in Assumption 1 (see in Appendix). However, in our paper, we directly estimate QM,D based on the concentration inequality for self-adjoint operators (Rudi & Rosasco, 2017; Lin & Cevher, 2020; 2018; Caponnetto & Yao, 2006), and prove that QM,D = O ( 1/(λ|D|) + √ 1/(λ|D|) ) (see in Proposition 6). Thus, our
estimation of QM,D is √ N (λ) tighter than that in (Li et al., 2019b). This is the one of the key reasons why we can substantially relax the restriction on the number of local machines compared with (Li et al., 2019b); (b) We not only present the bounds in expectation, but also in probability. To derive the tight bound in probability, we provide new decompositions of ‖f̄0M,D,λ − fM,D,λ‖ and ‖w̄0M,D,λ−wM,D,λ‖, please see Proposition 9 in detail. As far as we know, these decompositions are novel; (c) We also consider the Newton Raphson iteration-based communication strategy. To derive the improved high-probability bounds with communication, we introduce a novel decomposition of ‖f̄ tM,D,λ − fM,D,λ‖ρ (see in Proposition 10).
Compared with (Lin et al., 2020). (a) In (Lin et al., 2020), they also considered a communication strategy to enlarge the number of local machines. At first it seems that it only need to communicate the gradient information, but it should be noted that the gradient information is based on an operator representation (see Eq.(7) in (Lin et al., 2020)), which is usually infeasible in practice. The authors present a realization for the proposed strategy by communicating the data among each local machine, see in Appendix B (page 34 in (Lin et al., 2020), step 1). Thus, the data privacy of each local machine is difficult to be protected. Furthermore, since it requires communicating data Dj , j = 1, . . . ,m, among each local machine, for each iteration, the communication complexity of each local machine is O(|D|d), which is too high for large scale data sets. However, the communication strategy proposed in this paper only requires communicating the gradient gM,Dj ,λ(w̄ t−1 M,D,λ) and the model parameters βt−1j , rather than the data, therefore our proposed strategy do better on privacy protection. Moreover, the communication complexity is only O(M) for each local machine, M |D|, which is suitable for large scale data sets; (b) At first it seems that the proof techniques of (Lin et al., 2020) can be easily extended to our paper, but it is not true. If we use the same proof technical of (Lin et al., 2020), we can only obtain that ‖f M,D,λ − fM,λ‖ρ = ‖(LM,D + λI)−1(LM,D − LM )(fρ − fM,λ)‖ = O((1/( √ λ|D|) + √ N (λ)/|D|)‖fρ − fM,λ‖ρ/ √ λ), where f M,D,λ, fM,λ and fρ are defined in Definition 12. Combing with Proposition 2, one can only obtain that ‖fM,D,λ − fM,λ‖ρ = O (( 1/λ|D|+ √ N (λ)/λ|D| ) ‖fρ − fM,λ‖ρ ) . However, in our paper, we introduce new decompositions of ‖f M,D,λ − fM,λ‖ρ and ‖fM,D,λ − f M,D,λ‖ρ (see Propositions 2 and 3 in detail), and further obtain that ‖fM,D,λ − fM,λ‖ρ = O (( 1/( √ λ|D|) + √ N (λ)/|D| ) ‖fρ −
fM,λ‖ρ ) (see in Proposition 4), which is Ω(1/ √ λ) tighter than the directly use of the techniques of (Lin et al., 2020). The novel decompositions ‖f M,D,λ − fM,λ‖ρ and ‖fM,D,λ − f M,D,λ‖ρ are the key reasons why we can guarantee the optimal performance even under m = Ω( √ |D|). We only give an example here, but the novel decompositions have also been embedded in Propositions 5, 9, 10, etc.
Compared with (Rudi & Rosasco, 2017). (a) We study the statistical properties of the combination of distributed learning and random features, but in (Rudi & Rosasco, 2017), they only consider random features. As mentioned above, to estimate the tight bound of the distributed error, we introduce a novel decomposition of ‖f̄0M,D,λ − fM,D,λ‖ to derive the tight bound in expectation, and novel decompositions of ‖f̄0M,D,λ− fM,D,λ‖ and ‖f̄ p M,D,λ− fM,D,λ‖ to derive tight bounds in probability; (b) The combination not only brings distributed error, but also brings some other problems. If you compare the proofs of (Rudi & Rosasco, 2017) with these of our paper in detail, you can find that the decompositions of (Rudi & Rosasco, 2017) and ours are very different.
Overall, we improve the existing state-of-art bounds in expectation, and provide novel communicationbased distributed bounds with RF in probability. Moreover, we introduce some novel techniques and decompositions to substantially relax the restriction on the number of local machines, which are non-trivial extensions of (Li et al., 2019b; Lin et al., 2020; Rudi & Rosasco, 2017).

heading: 5 EXPERIMENTS
text: In this section, we validate our theoretical findings by performing experiments on both simulated and real datasets.
Numerical Experiments. Inspired by numerical experiments in (Rudi & Rosasco, 2017; Li et al., 2019e), we consider a spline kernel of order q: K2q(x,x′) = 1 + ∑∞ k=1 cos(2πk(x− x′))/(k2q). If
the marginal distribution of X is uniform on [0,1], then K2q(x,x′) = ∫ 1
0 ψ(x,ω)ψ(x,ω)%(ω)dω,
where ψ(x,ω) = Kq(x,ω) and %(ω) is also uniform on [0,1]. The random features of the spline kernel are φM (x) = (ψ(x,ω1), . . . , ψ(x,ωM ))T/ √ M. According to Theorem 1, 2 and 3, we set
the size of the random features to be M = √ |D|, and fine tune λ around |D|−1/2 using 5-fold cross validation1, the tuned set is {2−5, 2−3, . . . 25}|D|−1/2. We let the target function f∗ be a Gaussian random variable with mean µ = Kt(x, 0) and variance σ2 = 0.01.
We generate 10000 samples for training and 10000 samples for testing. We use the exact KRR as a baseline, which trains all samples in a batch. We compare our proposed DKRR-RF-CM (p = 2, 4, 8) with KRR and DKRR-RF. We repeat the training 5 times and estimate the averaged error on testing data. The mean square error on the test set with different partitions is given in Figure 1(a, b), which can be summarized as follows: 1) When m is not too large, the distributed methods (DRKK-RF and DRKK-RF-CM) are always comparable to original KRR. There exists an upper bound of m, when larger than it, the error increases dramatically and is far from the original KRR. This verifies the theoretical statement in Theorem 1, 2 and 3; 2) The upper bound m of DKRR-RF-CM is much larger than DKRR-RF. This result is aligned with Theorem 3, which demonstrates that the proposed communication strategy can enlarge the upper bound m; 3) The upper bound m of DKRR-RF-CM monotonically increases with the number of communications, which verifies Theorem 3.
Real Data. In this experiment, we consider the performance on real data. We use 6 publicly available datasets from LIBSVM Data2. The empirical evaluations with Gaussian kernel, exp(−‖x− x′‖2/σ), are given in Figure 2, where the optimal σ and λ are selected by 5-fold cross-validation, σ ∈ {2i, i = −10,−8, . . . , 10}, {2−5, 2−3, . . . 25}|D|−1/2, and the number of random features is 2 √ |D|.
1Because of the selection of the optimal λ using the 5-fold cross validation, the computational complexity should be enlarged. However, it should be noted that even for the plain methods, tuning the optimal λ is also required, which will enlarge the computational complexity as well.
2http://www.csie.ntu.edu.tw/∼cjlin/libsvm.
From Figure 2, one can find that: 1) our DKRR-RF-CM are better than the original DKRR-RF on all data-sets; 2) the larger the iterations of the communication, the better the performance. The above results demonstrate that our communication-based DKRR-RF is effective.

heading: 6 CONCLUSION
text: In this paper, we study the generalization properties of the combination of distributed learning and random features for ridge regression. We first improve the existing results of divide-and-conquer KRR with random features in expectation. Then, beyond the expectation, we derive generalization error bound in probability. Finally, we propose a novel effective communication strategy to further improve the learning performance of the combination method, and demonstrate the power of communications via both theoretical assessments and numerical experiments. Our results may open several venues for both theoretical and empirical work: (a) combine the approach with gradient algorithms such as preconditioned conjugate gradient Avron et al. (2017a) and multi-pass SGD (Carratino et al., 2018; Lin & Cevher, 2018; 2020); (b) replace synchronous distributed methods with asynchronous ones (Suresh et al., 2017); (c) consider the loss functions other than quadratic loss (Li et al., 2019e).

heading: Acknowledgment
text: This work was supported in part by the National Natural Science Foundation of China NO. 62076234 and the Beijing Outstanding Young Scientist Program NO. BJJWZYJH012019100020098.

heading: A APPENDIX: FAST LEARNING RATES
text: In the section, we will show that the fast learning rates can be achieved under favorable conditions.
Let
L2ρX =
{ f : X → R ∣∣∣ ∫ X f2(x)dρX ≤ ∞ } be the square integrable space, ‖ · ‖2ρ be its norm. Denote the integral operator (Smale & Zhou, 2007) LK by
LKf = ∫ X K(x, ·)f(x)dρX ,∀f ∈ L2ρX .
Assumption 1. For λ > 0, N (λ) is the effective dimension of the integral operator LK defined as N (λ) = Tr ( (LK + λI) −1LK ) , where Tr is the trace. Assuming there exists a constant c ≥ 1, such that
N (λ) ≤ cλ−γ , γ ∈ [0, 1]. (7)
The effective dimension is a common assumption within the framework of learning theory (Caponnetto & Vito, 2007; Smale & Zhou, 2005; Rudi & Rosasco, 2017), which is used to measure the complexity of the hypothesis space. It is always satisfied for γ = 1 and c = κ. Equation 7 can control the variance of the estimator and is equivalent to the classic entropy and covering number conditions. In particular, it holds if the eigenvalues of integral operator LK decay as i−1/γ , which is satisfied by the popular Gaussian and polynomial kernel functions. More details can be seen in (Steinwart & Christmann, 2008; Caponnetto & Vito, 2007; Rudi & Rosasco, 2017). Assumption 2. Let fρ(x) = ∫ Y ydρ(y|x) be the regression function, ρ(y|x) be the conditional distribution at x induced by ρ. For 12 ≤ r ≤ 1, assume there exists a g ∈ L 2 ρX such that
fρ(x) = L r Kg(x), (8)
where LrK is the rth power of LK .
Regression function fρ is the best function in L2ρX , which is the primary objective in regression problem. Assumption 2 is used to measure the complexity of the regression function fρ, which is commonly used in approximation theory (Caponnetto & Vito, 2007). Equation 8 can be used to control the bias of the estimator, it requires the expansion of the regression function fρ having coefficients that decay faster than the eigenvalues of integral operator LK . The larger the value of r, the faster the coefficients decay. The case r = 1/2 means that fρ ∈ HK . More detail can be seen in (Rudi & Rosasco, 2017; Smale & Zhou, 2007). Theorem 4. Suppose ψ is continuous, such that |ψ(x,ω)| ≤ τ almost surely, τ ∈ [1,∞) and |y| ≤ ζ . Under Assumptions 1-2 with r ∈ [1/2, 1], γ ∈ [0, 1], if λ = Ω(|D|− 1 2r+γ ), |D1| = . . . |Dm|,
m . |D| 2r+γ−1 2r+γ and M & |D| 1+γ(2r−1) 2r+γ , (9)
then, for every δ ∈ (0, 1], with probability at least 1− δ, we have
E [ ‖f̄0M,D,λ − fρ‖2ρ ] = O ( |D|− 2r 2r+γ log2 1
δ
) .
The bound above is the same as the original KRR estimator and is optimal in a minimax sense (Caponnetto & Vito, 2007; Lin et al., 2017; Chang et al., 2017b). In the best case, when r = 1 and γ = 0, the rate O(1/|D|) can be achieved by Ω( √ |D|) random features and Ω( √ |D|) partitions. In the worst case, that is r = 1/2 and γ = 1, which has been covered in Theorem 1. Theorem 5. Suppose ψ is continuous, such that |ψ(x,ω)| ≤ τ almost surely, τ ∈ [1,∞) and |y| ≤ ζ . Under Assumptions 1-2 with r ∈ [1/2, 1], γ ∈ [0, 1], if λ = Ω(|D|− 1 2r+γ ), |D1| = . . . |Dm|,
m . |D| 2r+γ−1 4r+2γ and M & |D| 1+γ(2r−1) 2r+γ , (10)
then, for every δ ∈ (0, 1], with probability at least 1− δ, we have ‖f̄0M,D,λ − fρ‖2ρ = O ( |D|− 2r 2r+γ log2 1
δ
) .
One can see that the upper bound of m is |D| 2r+γ−1 4r+2γ , which is stricter than that of Theorem 4. In the best case, when r = 1 and γ = 0, the rate O(1/|D|) can be achieved considering Ω( √ |D|) random
features and Ω(|D|1/4) partitions. For the basic setting, that is r = 1/2 and γ = 1, which has been given in Theorem 2. Theorem 6. Suppose ψ is continuous, such that |ψ(x,ω)| ≤ τ almost surely, τ ∈ [1,∞) and |y| ≤ ζ . Under Assumptions 1-2 with r ∈ [1/2, 1], γ ∈ [0, 1], if λ = Ω(|D|− 1 2r+γ ), |D1| = . . . |Dm|,
m . |D| (2r+γ−1)(p+1) (2r+γ)(p+2) and M & |D| 1+γ(2r−1) 2r+γ log 1
δ , (11)
then, for any δ ∈ (0, 1], with probability at least 1− δ,
‖f̄pM,D,λ − fρ‖ 2 ρ = O
( |D|− 2r 2r+γ logp+2 1
δ
) .
One can see that the communication can relax the restriction on the number of partitions. As p→∞, the partitions can reach Ω(|D| 2r+γ−1 2r+γ ), which is the same as Theorem 4.
Remark 5. In this paper, we focus on enlarging the number of local machines. In (Rudi & Rosasco, 2017; Rudi et al., 2018), they have proved that when generating random features in a data-dependent manner, fewer random features are required to obtain optimal learning. Thus, if we adopt the data-dependent manner to generate random features, we can further improve the performance of DKRR-RF-CM with fewer random features.

heading: B APPENDIX: NOTATION AND PRELIMINARY
text: In this paper we denote the operator norm by ‖ · ‖, and the square integrable norm by ‖ · ‖ρ. Definition 1.
SM : RM → L2ρX , (SMw)(x) = 〈w,φM (x)〉,
S∗M : L 2 ρX → R M , S∗Mg = ∫ X φM (x)g(x)dρX (x),
S∗M,D : L 2 ρX → R M , S∗M,Dg = 1 |D| ∑
xj∈DX
φM (xj)g(xj),
CM : RM → RM , CM = ∫ X φM (x)φM (x) TdρX (x),
CM,D : RM → RM , CM,D = 1 |D| ∑
xj∈DX
φM (xj)φM (xj) T.
Lemma 1. CM and CM,D are self-adjoint and positive operators, with spectrum is [0, τ2]. Moreover we have CM = S∗MSM and CM,D = ΦM,DΦ T M,D = S ∗ M,DSM .
Proof. CM and CM,D are self-adjoint and positive operators, with spectrum is [0, τ2], and
CM = S ∗ MSM , CM,D = ΦM,DΦ T M,D
can be directly obtained from Caponnetto & Vito (2007); Smale & Zhou (2005; 2007); Rosasco et al. (2010); Rudi & Rosasco (2017); Lin et al. (2020).
In the following, we prove that CM,D = S∗M,DSM . From the definitions of SM , S ∗ M,D and CM,D in Definition 1, we have
∀β ∈ RM , (SMβ)(·) = 〈β,φM (·)〉 = φM (·)Tβ, and thus we can obtain that,
S∗M,DSMβ = 1 |D| ∑
xj∈DX
φM (xj)φM (xj) Tβ = CM,Dβ.
So, the Equation CM,D = S∗M,DSM holds.
Definition 2.
fM,D,λ = w T M,D,λφM (·),wM,D,λ = arg min
w∈RM
{ 1
|D| ∑ zi∈D ( wTφM (xi)− yi )2 + λ‖w‖2
} ;
f M,D,λ = w T M,D,λφM (·),w M,D,λ = arg min
w∈RM
{ 1
|D| ∑ zi∈D ( wTφM (xi)− fρ(xi) )2 + λ‖w‖2
} ;
fM,λ = w T M,λφM (·),wM,λ = arg min
w∈RM {∫ X ( wTφM (x)− fρ(x) )2 dρX (x) + λ‖w‖2 } ;
fλ = arg min f∈HK {∫ X (f(x)− fρ(x))2 dρX (x) + λ‖f‖2K } .
(12) Remark 6 (From Caponnetto & Vito (2007); Smale & Zhou (2007); Li et al. (2019b); Lin et al. (2020)).
fM,D,λ = SMwM,D,λ,wM,D,λ = (CM,D + λI) −1ΦM,DȳD;
f M,D,λ = SMw M,D,λ,w M,D,λ = (CM,D + λI) −1S∗M,Dfρ;
fM,λ = SMwM,λ,wM,λ = (CM + λI) −1S∗Mfρ.
(13)
Definition 3. The maximum random feature dimension is denoted as N∞(λ) = sup
ω∈Ω ‖(LK + λI)−1/2ψ(·,ω)‖2ρ, λ > 0.
Remark 7 (From Rudi & Rosasco (2017)). N∞(λ) ≤ τ2λ−1 is always satisfied for every λ > 0.

heading: C APPENDIX: PROOF OF THEOREM 4

heading: C.1 APPENDIX: ERROR DECOMPOSITION FOR DKRR-RF IN EXPECTATION
text: Proposition 1. Let f̄0M,D,λ, and f M,D,λ, fM,λ, fλ be defined by 3 and 12, respectively. Then, we have E [ ‖f̄0M,D,λ − fρ‖2ρ ] ≤ 3 [ ‖fM,λ − fλ‖2ρ ] + 3 [ ‖fλ − fρ‖2ρ ] + 3
m∑ j=1 |Dj |2 |D|2 E [ ‖fM,Dj ,λ − fM,λ‖2ρ ] + 3 m∑ j=1 |Dj | |D| E [ ‖f M,Dj ,λ − fM,λ‖ 2 ρ ] .
Proof. Let f̄0M,D,λ, and f M,D,λ, fM,λ be defined by 3 and 12, respectively. According to the Proposition 5 of (Chang et al., 2017b) or Lemma 20 of (Li et al., 2019b), we have
E [ ‖f̄0M,D,λ − fM,λ‖2ρ ] ≤
m∑ j=1 |Dj |2 |D|2 E [ ‖fM,Dj ,λ − fM,λ‖2ρ ] + m∑ j=1 |Dj | |D| E [ ‖f M,D,λ − fM,λ‖2ρ ] .
(14)
Note that (a+ b+ c)2 ≤ 3a2 + 3b2 + 3c2, ∀a, b, c ≥ 0. Thus, we have E [ ‖f̄0M,D,λ − fρ‖2ρ ] = E [ ‖f̄0M,D,λ − fM,λ + fM,λ − fλ + fλ − fρ‖2ρ ] ≤ 3E [ ‖f̄0M,D,λ − fM,λ‖2ρ ] + 3 [ ‖fM,λ − fλ‖2ρ ] + 3 [ ‖fλ − fρ‖2ρ ] .
Combining the above inequality and Eq. 14, we can prove this proposition.
Proposition 2. The follows hold: √ λ‖wM,D,λ −w M,D,λ‖ ≤ JM,D(RM,D +KM,D)
and ‖fM,D,λ − f M,D,λ‖ρ ≤ J 2M,D(RM,D +KM,D).
where RM,D := ∥∥(CM + λI)−1/2(ΦM,DȳD − S∗Mfρ)∥∥, JM,D := ‖(CM,D + λI)−1/2(CM + λI)1/2‖, KM,D := ‖(CM + λI)−1/2(S∗Mfρ − S∗M,Dfρ)‖.
Proof. From 13, we know that wM,D,λ = (CM,D + λI)−1ΦM,DȳD and w M,D,λ = (CM,D + λI)−1S∗M,Dfρ, so we have
wM,D,λ −w M,D,λ = (CM,D + λI)−1(ΦM,DȳD − S∗M,Dfρ)
=(CM,D + λI) −1/2(CM,D + λI) −1/2(CM + λI) 1/2(CM + λI) −1/2(ΦM,DȳD − S∗M,Dfρ). (15) Note (CM,D + λI)−1/2 is a self-adjoint and positive operator, so ‖(CM,D + λI)−1/2‖ ≤ 1/ √ λ, thus, we can obtain that
‖wM,D,λ −w M,D,λ‖ ≤ 1√ λ JM,D‖(CM + λI)−1/2(ΦM,DȳD − S∗M,Dfρ)‖
= 1√ λ JM,D‖(CM + λI)−1/2(ΦM,DȳD − S∗M,Dfρ + S∗M,Dfρ − S∗M,Dfρ)
≤ 1√ λ JM,D (RM,D +KM,D) .
(16)
Note that fM,D,λ − f M,D,λ = SM (wM,D,λ −w M,D,λ), by 15, we have fM,D,λ − f M,D,λ = SM (wM,D,λ −w M,D,λ)
=SM (CM + λI) −1/2(CM + λI) 1/2(CM,D + λI) −1/2(CM,D + λI) −1/2(CM + λI) 1/2
(CM + λI) −1/2(ΦM,DȳD − S∗M,Dfρ + S∗M,Dfρ − S∗M,Dfρ).
(17)
Note that ‖SM (CM + λI)−1/2‖ = ‖(CM + λI)−1/2S∗MSM (CM + λI)−1/2‖1/2 = ‖(CM + λI)−1/2CM (CM + λI)−1/2‖1/2 ≤ 1. So, by Eq 17, we have
‖fM,D,λ − f M,D,λ‖ρ ≤ J 2M,D (RM,D +KM,D) .
Proposition 3. The follows hold: √ λ‖w M,D,λ −wM,λ‖ ≤ ‖fM,λ − fρ‖ρ + JM,D‖fM,λ − fρ‖ρ
and ‖f M,D,λ − fM,λ‖ρ ≤ (JM,D + J 2M,D)‖fM,λ − fρ‖ρ,
where JM,D := ‖(CM,D + λI)−1/2(CM + λI)1/2‖.
Proof. By Remark 6, we have w M,D,λ −wM,λ = (CM,D + λI)−1S∗M,Dfρ − (CM + λI)−1S∗Mfρ
=(CM,D + λI) −1[S∗M,Dfρ − S∗Mfρ] + [(CM,D + λI)−1 − (CM + λI)−1]S∗Mfρ.
Note that for any self-adjoint and positive operators A and B, A−1 −B−1 = A−1(B −A)B−1, A−1 −B−1 = B−1(B −A)A−1, (18)
so we have w M,D,λ −wM,λ
=(CM,D + λI) −1[S∗M,Dfρ − S∗Mfρ] + (CM,D + λI)−1(CM − CM,D)wM,λ
From Lemma 1, we know that CM = S∗MSM and CM,D = ΦM,DΦ T M,D = S ∗ M,DSM , thus we can obtain that w M,D,λ −wM,λ
=(CM,D + λI) −1[S∗M,Dfρ − S∗Mfρ] + (CM,D + λI)−1(S∗MSMwM,λ − S∗M,DSMwM,λ) =(CM,D + λI) −1[S∗M,Dfρ − S∗M,DSMwM,λ] + (CM,D + λI)−1[S∗MSMwM,λ − S∗Mfρ] =(CM,D + λI) −1[S∗M,Dfρ − S∗M,DfM,λ] + (CM,D + λI)−1[S∗MfM,λ − S∗Mfρ] =(CM,D + λI) −1S∗M,D[fρ − fM,λ] + (CM,D + λI)−1S∗M [fM,λ − fρ]. (19)
Thus, we have ‖w M,D,λ −wM,λ‖ ≤ ( ‖(CM,D + λI)−1S∗M,D‖+ ‖(CM,D + λI)−1S∗M‖ ) ‖fM,λ − fρ‖ρ. (20)
Note that
‖(CM,D + λI)−1/2S∗M,D‖ ≤ ‖(CM,D + λI)−1/2CM,D(CM,D + λI)−1/2‖1/2 ≤ 1
and
‖(CM,D + λI)−1/2S∗M‖ =‖(CM,D + λI)−1/2(CM + λI)1/2(CM + λI)−1/2S∗M‖ ≤JM,D‖(CM + λI)−1/2S∗M‖ ≤ JM,D,
since ‖(CM + λI)−1/2S∗M‖ = ‖(CM + λI)−1/2CM (CM + λI)−1/2‖1/2 ≤ 1. Substituting the above two inequalities into Eq. 20, we have
‖w M,D,λ −wM,λ‖ ≤ 1√ λ (1 + JM,D)‖fM,λ − fρ‖ρ,
which prove the first result of Proposition 3.
In the following, we will prove the second result. By Eq. 19, we have
f M,D,λ − fM,λ = SM (w M,D,λ −wM,λ) =SM (CM,D + λI) −1S∗M,D[fρ − fM,λ] + SM (CM,D + λI)−1S∗M [fM,λ − fρ]
=SM (CM + λI) −1/2(CM + λI) 1/2(CM,D + λI) −1/2(CM,D + λI) −1/2S∗M,D[fρ − fM,λ]
+ SM (CM + λI) −1/2(CM + λI) 1/2(CM,D + λI) −1/2(CM,D + λI) −1/2(CM + λI) 1/2
(CM + λI) −1/2S∗M [fM,λ − fρ]
Note that ‖SM (CM + λI)−1/2‖ = ‖(CM + λI)−1/2CM (CM + λI)−1/2‖1/2 ≤ 1, ‖(CM,D + λI)−1/2S∗M,D‖ = ‖(CM,D + λI)−1/2CM,D(CM,D + λI)−1/2‖1/2 ≤ 1, and ‖(CM + λI)−1/2S∗M‖ = ‖(CM + λI)−1/2CM (CM + λI)−1/2‖1/2 ≤ 1, so we have
‖f M,D,λ − fM,λ‖ρ ≤ (JM,D + J 2M,D)‖fM,λ − fρ‖ρ. (21)
Proposition 4. The follows hold: √ λ‖wM,D,λ −wM,λ‖ ≤ JM,D(RM,D +KM,D) + (1 + JM,D)‖fM,λ − fρ‖ρ,
and
‖fM,D,λ − fM,λ‖ρ ≤ J 2M,D(RM,D +KM,D) + (JM,D + J 2M,D)‖fM,λ − fρ‖ρ.
where RM,D := ∥∥(CM + λI)−1/2(ΦM,DȳD − S∗Mfρ)∥∥, JM,D := ‖(CM,D + λI)−1/2(CM + λI)1/2‖, KM,D := ‖(CM + λI)−1/2(S∗Mfρ − S∗M,Dfρ)‖.
Proof. Note that
‖wM,D,λ −wM,λ‖ ≤ ‖wM,D,λ −w M,D,λ‖+ ‖w M,D,λ −wM,λ‖
and
‖fM,D,λ − fM,λ‖ρ ≤ ‖fM,D,λ − f M,D,λ‖ρ + ‖f M,D,λ − fM,λ‖ρ.
Combining Propositions 2 and 3, we can prove this result.

heading: C.2 PROOF OF THEOREM 4
text: Lemma 2 (Lemma 23 in (Li et al., 2019b), can be also seen in (Rudi & Rosasco, 2017)). For δ ∈ (0, 1] and λ > 0, when
M ≥ Ω (( N (λ) λ )2r−1( N∞(λ) log 1 λ )2−2r ∨ (N∞(λ)) log 1 λδ ) .
Then, with probability at least 1− δ, we have
‖fM,λ − fλ‖2ρ ≤ cλ2r,
where c is a constant.
Lemma 3 (Theorem 4 in (Smale & Zhou, 2005)). Under Assumption 2, for r ∈ [1/2, 1], we have
‖fλ − fρ‖2ρ ≤ cλ2r,
where c is a constant.
Lemma 4 (Lemma 6 in (Rudi & Rosasco, 2017)). For δ ∈ (0, 1], with probability at least 1− δ, we have
RM,D := ∥∥∥(CM + λI)−1/2 (ΦM,DȳD − S∗Mfρ)∥∥∥ = O (( 1√ λ|D| + √ NM (λ) |D| ) log 1 δ ) ,
where NM (λ) := Tr ( (LM + λI) −1LM ) , LM is the integral operator associated with the approxi-
mate kernel function KM , (LMf)(x) = ∫ X KM (x,x ′)f(x′)dρX (x ′).
Lemma 5 (Proposition 10 in (Rudi & Rosasco, 2017)). For any δ ∈ (0, 1], M ≥ Ω ( N∞(λ) log 1λδ ) , then with probability at least 1− δ,
|NM (λ)−N (λ)| ≤ 1.55N (λ), (22) where NM (λ) := Tr ( (LM + λI) −1LM ) .
Lemma 6 (Lemma E.2, (Blanchard & Krämer, 2010)). For any self-adjoint and positive semidefinite operators A and B, if there exists η > 0 such that the following inequality holds
‖(A+ λI)−1/2(B −A)(A+ λI)−1/2‖ ≤ 1− η,
then
‖(A+ λI)1/2(B + λI)−1/2‖ ≤ 1√ η .
Proposition 5. For δ ∈ (0, 1], with probability at least 1− δ, we have
KM,D := ‖(CM + λI)−1/2(S∗Mfρ − S∗M,Dfρ)‖ ≤ 2τζ log 1δ 3|D| √ λ + 2ζ √ NM (λ) |D| ,
where NM (λ) := Tr ( (LM + λI) −1LM ) .
The proof closely follows the proof of Lemma 6 in (Rudi & Rosasco, 2017).
Proof. Define µi = (CM + λI)1/2S∗Mfρ − (CM + λI)1/2φM (xi)fρ(xi). Note that
(CM + λI) 1/2(S∗Mfρ − S∗M,Dfρ) =
1
|D| |D|∑ i=1 µi.
Since µ1, . . . , µ|D| are independent and identically distributed random vector, and Eµi = ∫ X (CM + λI) 1/2φM (x)fρ(x)dρX − ∫ X (CM + λI) 1/2φM (xi)fρ(xi)dρX = 0.
To apply the Bernstein inequality (Arcones, 1995; Rudi & Rosasco, 2017) for random vectors, we need to bound their moments. Note that
‖(CM + λI)1/2φM (x)fρ(x)‖ ≤ τζ√ λ ,
and
Eµ2i ≤ 2 ∫ X ‖(CM + λI)1/2φM (x)‖2‖fρ(x)‖2ρdρX (x)
≤ 2ζ2 ∫ X ‖(CM + λI)1/2φM (x)‖2dρX (x) ≤ 2ζ2NM (λ).
Thus, using the Bernstein inequality (Arcones, 1995; Rudi & Rosasco, 2017), for any δ ∈ (0, 1], with 1− δ, we have
‖(CM + λI)−1/2(S∗Mfρ − S∗M,Dfρ)‖ ≤ 2τζ log 1δ 3|D| √ λ + 2ζ √ NM (λ) |D| .
Lemma 7 (Lemma 2 in (Smale & Zhou, 2007)). Let H be a Hilbert space and ξ be a random variable on Z, ρ with values inH. Assume ‖ξ‖ ≤M <∞ almost surely. Denote σ(ξ) = E(ξ2). Let {zi}ni=1 be independent random drawers of ρ. For any 0 < δ < 1, with confidence 1− δ,∥∥∥∥∥ 1n m∑ i=1 [ξ − E(ξ)] ∥∥∥∥∥ ≤ 2M log(2/δ)n + √ 2‖σ(ξ)‖2 log(2/δ) n .
Proposition 6. For any δ > 0, with probability at least 1− δ, we have
QM,D := ‖(CM + λI)−1/2(CM − CM,D)(CM + λ)−1/2‖
≤ 2 log 2(2/δ)(N∞(λ) + 1)
|D| +
√ 2 log(2/δ)(N∞(λ) + 1)
|D| .
and
LM,D := ‖(CM + λI)−1(CM − CM,D)‖
≤ 2 log 2(2/δ)(N∞(λ) + 1)
|D| +
√ 2 log(2/δ)(N∞(λ) + 1)
|D| .
where N∞(λ) = supω∈Ω ‖(LK + λI)−1/2ψ(·,ω)‖2ρ, c1 and c2 are two constants.
To prove Proposition 6, we first prove the following lemma (a similar technique can be found in (Hsu et al., 2012; Rudi et al., 2013; Caponnetto & Yao, 2006; Rudi & Rosasco, 2017)).
The similar result for the matrix case was first proved in (Hsu et al., 2012), and later was extended to the operator case in (Rudi et al., 2013; Rudi & Rosasco, 2017).
Lemma 8. Let ζ1, . . . , ζn with n ≥ 1, be i.i.d random vectors on a separable Hilbert spaces H such that H = Eζ ⊗ ζ is a trace class, and for any λ there exists N∞(λ) < ∞ such that 〈ζ, (H + λI)−1ζ〉 ≤ N∞(λ). Denote Hn as 1n ∑n i=1 ζi ⊗ ζi. Then for any δ ≥ 0, with probability at least 1− 2δ, the following holds∥∥∥(H + λI)−1/2(H −Hn)(H + λ)−1/2∥∥∥ ≤ 2 log2(2/δ)(N∞(λ) + 1) n + √ 2 log(2/δ)(N∞(λ) + 1) n .
Proof. Let Hλ = H + λI , η = H −1/2 λ HH −1/2 λ , ξi = η −H −1/2 λ ζi ⊗H −1/2 λ ζi. One can see that Eξi = 0. Note that
‖η −H−1/2λ ζi ⊗H −1/2 λ ζi‖ ≤ ‖η‖+ 〈ζi, H −1/2 λ ζi〉 ≤ 1 +N∞(λ),
and
‖E[ξ2i ]‖ = ∥∥∥E [〈ζi, H−1/2λ ζi〉H−1/2λ ζi ⊗H−1/2λ ζi]−H−2λ H2∥∥∥
≤ N∞(λ) ∥∥∥E[H−1/2λ ζi ⊗H−1/2λ ζi]∥∥∥+ ∥∥H−2λ H2∥∥
≤ N∞(λ)‖H−1/2λ H‖+ ∥∥H−2λ H2∥∥ ≤ N∞(λ) + 1.
Thus, substituting the above two inequalities to Lemma 7 (Lemma 2 in (Smale & Zhou, 2007)), which finishes the proof.
Proof of Proposition 6. Since CM is self-adjoint operator, so we have
‖(CM + λI)−1(CM − CM,D)‖ = ‖(CM + λI)−1/2(CM − CM,D)(CM + λI)−1/2‖.
According to Lemma 8 with υi = φM (xi), we can obtain this result.
Proposition 7. If |D| ≥ 32 log(2/δ)(1+N∞(λ)), then for any δ > 0, with probability at least 1−δ, we have
JM,D := ‖(CM,D + λI)−1/2(CM + λI)1/2‖ ≤ √ 2.
Proof. From Proposition 6, we know that if |D| ≥ 32 log(2/δ)(1 +N∞(λ)), then
‖(CM + λI)−1/2(CM,D − CM )(CM + λ)−1/2‖ ≤ 1
2 .
Combining the above inequality and Lemma 6, we can prove this result.
Proposition 8. If δ ∈ (0, 1], and |D| ≥ Ω(N∞(λ)), then with 1− δ, we have ‖fM,D,λ − fM,λ‖ρ = O ( ΥM,D,λ log 1
δ + ‖fM,λ − fλ‖ρ + ‖fλ − fρ‖ρ
) ,
where ΥM,D,λ := 1√λ|D| + √ N (λ) |D| .
Proof. From Proposition 4, we have
‖fM,D,λ − fM,λ‖ρ ≤ J 2M,D(RM,D +KM,D) + (JM,D + J 2M,D)‖fM,λ − fρ‖ρ.
Thus, from Lemmas 4, 5, and Propositions 5, 7, we know that if |D| ≥ Ω(N∞(λ)), we can prove this result.
Proof of Theorem 4. According to Proposition 1, we have E [ ‖f̄0M,D,λ − fρ‖2ρ ] ≤ 3 [ ‖fM,λ − fλ‖2ρ ] + 3 [ ‖fλ − fρ‖2ρ ] + 3
m∑ j=1 |Dj |2 |D|2 E [ ‖fM,Dj ,λ − fM,λ‖2ρ ] + 3 m∑ j=1 |Dj | |D| E [ ‖f M,Dj ,λ − fM,λ‖ 2 ρ ] .
Substituting Lemmas 2,3, Proposition 3, 7, 8 into the above inequality, one can see that if M ≥ Ω (( N 2r−1(λ) λ2r−1 ) N∞(λ)2−2r ∨N∞(λ) ) , and |Dj | ≥ Ω(N∞(λ)),
with confidence 1− δ, we have
E [ ‖f̄0M,D,λ − fρ‖2ρ ] = O λ2r + m∑ j=1 |Dj |2 |D|2 ( Υ2M,Dj ,λ log 2 1 δ ) + m∑ j=1 |Dj | |D| λ2r  , (23)
where ΥM,Dj ,λ = 1√ λ|Dj | + √ N (λ) |Dj | . If setting |D1| = . . . = |Dm|, λ = Ω(|D| − 12r+γ ), we have
ΥM,Dj ,λ = O (√ NM (λ) |Dj | + 1 |Dj | √ λ ) = O (√ m|D|− r 2r+γ +m|D|− 4r+2γ−1 4r+2γ ) . (24)
Note that if m ≤ Ω ( |D| 2r+γ−1 2r+γ ) and |D1| = . . . = |Dm|, and λ = Ω(|D|− 1 2r+γ ), we have
|Dj | = |D| m ≥ Ω
( |D| 1 2r+γ ) = Ω(N∞(λ)).
Thus, substituting 24 into 23, one can see that if m ≤ Ω ( |D| 2r+γ−1 2r+γ ) and
M ≥ Ω (( N 2r−1(λ) λ2r−1 ) N∞(λ)2−2r ∨N∞(λ) ) = Ω ( |D| 1+(2r−1)γ 2r+γ ) ,
with probability at least 1− δ, we have
E [∥∥f̄0M,D,λ − fρ∥∥2ρ] = O(|D|− 2r2r+γ log2 1δ + |D|−1 log2 1δ + |D|− 2r2r+γ ) = O ( |D|− 2r 2r+γ log2 1
δ
) .

heading: D APPENDIX: PROOF OF THEOREM 5

heading: D.1 APPENDIX: ERROR DECOMPOSITION FOR DKRR-RF IN PROBABILITY
text: Proposition 9. The follows hold:
‖w̄0M,D,λ −wM,D,λ‖ ≤ m∑ j=1 |Dj | |D| J 2M,D(QM,D +QM,Dj )‖wM,Dj ,λ −wM,λ‖ and
‖f̄0M,D,λ − fM,D,λ‖ ≤ m∑ j=1 |Dj | |D| J 2M,D ( QM,D +QM,Dj ) ( ‖fM,Dj ,λ − fM,λ‖ρ + √ λ‖wM,Dj ,λ −wM,λ‖ ) ,
where JM,D := ‖(CM,D + λI)−1/2(CM + λI)1/2‖ and QM,D := ‖(CM + λI)−1/2(CM − CM,D)(CM + λ) −1/2‖.
Proof. Note that wM,D,λ = (CM,D + λI)−1ΦM,DȳD, thus we have
w̄0M,D,λ −wM,D,λ
= m∑ j=1 |Dj | |D| (CM,Dj + λI) −1ΦM,Dj ȳDj − (CM,D + λI)−1ΦM,DȳD
= m∑ j=1 |Dj | |D| ( (CM,Dj + λI) −1 − (CM,D + λI)−1 ) ΦM,Dj ȳDj
= m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM,D − CM,Dj) (CM,Dj + λI)−1ΦM,Dj ȳDj
= m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM,D − CM,Dj)wM,Dj ,λ
= m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM,D − CM ) wM,Dj ,λ
+ m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM − CM,Dj)wM,Dj ,λ
= m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM,D − CM ) (wM,Dj ,λ −wM,λ)
+ m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM,D − CM ) wM,λ
+ m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM − CM,Dj)wM,Dj ,λ
= m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM,D − CM ) (wM,Dj ,λ −wM,λ)
+ m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM − CM,Dj) (wM,Dj ,λ −wM,λ).
(25)
Note that
m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM,D − CM ) (wM,Dj ,λ −wM,λ)
= m∑ j=1 |Dj | |D| (CM,D + λI) −1(CM + λI)(CM + λI) −1 (CM,D − CM ) (wM,Dj ,λ −wM,λ)
and m∑ j=1 |Dj | |D| (CM,D + λI) −1 (CM − CM,Dj) (wM,Dj ,λ −wM,λ)
= m∑ j=1 |Dj | |D| (CM,D + λI) −1(CM + λI)(CM + λI) −1 (CM − CM,Dj) (wM,Dj ,λ −wM,λ). Substituting the above equations into Eq. 25, we have
‖w̄0M,D,λ −wM,D,λ‖ ≤ m∑ j=1 |Dj | |D| J 2M,D(LM,D + LM,Dj )‖wM,Dj ,λ −wM,λ‖,
where LM,D := ‖(CM +λI)−1(CM −CM,D)‖. From Proposition 6, we know that LM,D = QM,D, so we have
‖w̄0M,D,λ −wM,D,λ‖ ≤ m∑ j=1 |Dj | |D| J 2M,D(QM,D +QM,Dj )‖wM,Dj ,λ −wM,λ‖,
which prove the first result of this proposition.
In the following, we will prove the second result of this proposition. Note that SM (w̄0M,D,λ − wM,D,λ) = f̄ 0 M,D,λ − fM,D,λ. According to 25, we have
f̄0M,D,λ − fM,D,λ = m∑ j=1 |Dj | |D| SM (CM,D + λI) −1 (CM,D − CM ) (wM,Dj ,λ −wM,λ)
+ m∑ j=1 |Dj | |D| SM (CM,D + λI) −1 (CM − CM,Dj) (wM,Dj ,λ −wM,λ)
:= m∑ j=1 |Dj | |D| (~1j + ~2j ).
(26)
Note that ~1j =SM (CM + λI)−1/2(CM + λI)1/2(CM,D + λI)−1/2(CM,D + λI)−1/2(CM + λI)1/2
(CM + λI) −1/2 (CM,D − CM ) (CM + λI)−1/2
(CM + λI) −1/2(CM + λI)(wM,Dj ,λ −wM,λ),
(27)
thus we have ‖~1j‖ρ ≤ J 2M,DQM,D‖SM (CM + λI)−1/2‖‖(CM + λI)−1/2(CM + λI)(wM,Dj ,λ −wM,λ)‖
≤ J 2M,DQM,D‖(CM + λI)−1/2(CM + λI)(wM,Dj ,λ −wM,λ)‖. (28)
Since ‖SM (CM + λI)−1/2‖ = ‖(CM + λI)−1/2CM (CM + λI)−1/2‖1/2 ≤ 1. So, we have ‖~1j‖ρ ≤ J 2M,DQM,D‖(CM + λI)−1/2(CM + λI)(wM,Dj ,λ −wM,λ)‖
= J 2M,DQM,D‖(CM + λI)−1/2(S∗MSM + λI)(wM,Dj ,λ −wM,λ)‖
≤ J 2M,DQM,D‖(CM + λI)−1/2S∗MSM (wM,Dj ,λ −wM,λ)‖
+ λJ 2M,DQM,D‖(CM + λI)−1/2(wM,Dj ,λ −wM,λ)‖
≤ J 2M,DQM,D‖(CM + λI)−1/2S∗M‖‖fM,Dj ,λ − fM,D,λ‖ρ + √ λJ 2M,DQM,D‖(wM,Dj ,λ −wM,λ)‖ ≤ J 2M,DQM,D(‖fM,Dj ,λ − fM,D,λ‖ρ + √ λ‖wM,Dj ,λ −wM,λ‖),
the last inequality uses the fact that
‖(CM + λI)−1/2S∗M‖ = ‖(CM + λI)−1/2S∗MSM (CM + λI)−1/2‖1/2 ≤ 1.
Similar as the above process, we can also obtain that
‖~2j‖ρ ≤ J 2M,DQM,Dj (‖fM,Dj ,λ − fM,D,λ‖ρ − √ λ‖wM,Dj ,λ −wM,λ‖).
Thus, using 26, we can prove this result.

heading: D.2 APPENDIX: PROOF OF THEOREM 5
text: Proof of Theorem 5. Combining Proposition 9 and Proposition 4, we have
‖f̄0M,D,λ − fM,D,λ‖ρ ≤ m∑ j=1 |Dj | |D| J 2M,D ( QM,D +QM,Dj ) (
(JM,Dj + J 2M,Dj )(RM,Dj +KM,Dj ) + (2JM,Dj + J 2 M,Dj + 1)‖fM,λ − fρ‖ρ
) .
(29)
From Propositions 7 and 8, one can see that if |Dj | ≥ Ω(N∞(λ)), and λ ≤ ‖LK‖, then for any δ > 0, with probability at least 1− δ,
‖f̄0M,D,λ − fM,D,λ‖ρ
=O  m∑ j=1 |Dj | |D| ( QM,D +QM,Dj )( ΥM,Dj ,λ log 1 δ + ‖fM,λ − fλ‖ρ + ‖fλ − fρ‖ρ ) , where ΥM,Di,λ =
1√ λ|Di| + √ N (λ) |Di| . Note that QM,D ≤ QM,Dj , and by Lemmas 2 and 3, so we
have
‖f̄0M,D,λ − fM,D,λ‖ρ = O  m∑ j=1 |Dj | |D| QM,DjΥM,Dj ,λ log 1 δ + λrQM,Dj  . (30) Note that
‖f̄0M,D,λ − fρ‖ρ = ‖f0M,D,λ − fM,D,λ + fM,D,λ − fM,λ + fM,λ − fλ + fλ − fρ‖ρ ≤‖f̄0M,D,λ − fM,D,λ‖ρ + ‖fM,D,λ − fM,λ‖ρ + ‖fM,λ − fλ‖ρ + ‖fλ − f‖ρ.
(31)
Combining E.q 30, Lemmas 2, 3, Proposition 8 and E.q 31, one can see that if M ≥ Ω (( N 2r−1(λ) λ2r−1 ) N∞(λ)2−2r ∨N∞(λ) ) , with probability 1− δ, we have
‖f̄0M,D,λ − fρ‖ρ = O  m∑ j=1 |Dj | |D| QM,DjΥM,Dj ,λ log 1 δ + ΥM,D,λ log 1 δ + λrQM,Dj + λr  . (32)
If |D1| = . . . = |Dm|, λ = Ω(|D|− 1 2r+γ ), we have ΥM,D,λ = O ( |D|− r 2r+γ + |D|− 4r+2γ−1 4r+2γ ) = O ( |D|− r 2r+γ ) ,
ΥM,Dj ,λ = O (√ m|D|− r 2r+γ +m|D|− 4r+2γ−1 4r+2γ ) ,
QM,Dj = O ( m|D|− 2r+γ−1 2r+γ + √ m|D|− 2r+γ−1 4r+2γ ) .
(33)
Thus, when m ≤ Ω ( |D| 2r+γ−1 4r+2γ ) , we have
ΥM,Dj ,λQM,Dj = O ( |D|− r 2r+γ ) ,
QM,Djλr = O ( |D|− r 2r+γ |D|− 2r−1+γ 8r+4γ ) = O ( |D|− r 2r+γ ) ,
|Dj | = |D| m ≥ Ω
( |D| 2r+γ+1 4r+2γ ) ≥ Ω ( |D| 1 2r+γ ) = Ω(N∞(λ)).
Thus, we have ‖f̄0M,D,λ − fρ‖2ρ = O ( |D|− 2r 2r+γ log2 1δ ) , which prove this result.

heading: E APPENDIX: PROOF OF THEOREM 6

heading: E.1 APPENDIX: ERROR DECOMPOSITION FOR DKRR-RF-CM IN PROBABILITY
text: Proposition 10. ‖f̄ tM,D,λ − fM,D,λ‖ρ
≤ ( m∑ i=1 |Dj | |D| (2J 2M,DjQM,D + 2J 2 M,DjQM,Dj )t ( ‖f̄0M,D,λ − fM,D,λ‖ρ + √ λ‖w̄0M,D,λ −wM,D,λ‖ ) ,
where JM,D := ‖(CM,D + λI)−1/2(CM + λI)1/2‖ and QM,D := ‖(CM + λI)−1/2(CM − CM,D)(CM + λ) −1/2‖.
Proof. Note that
wM,D,λ = w̄ t−1 M,D,λ − (CM,D + λI) −1 [ (CM,D + λI)w̄ t−1 M,D,λ −ΦM,DȳD ] ,
w̄tM,D,λ = w̄ t−1 M,D,λ − m∑ j=1 |Dj | |D| (CM,Dj + λI) −1 [ (CM,D + λI)w̄ t−1 M,D,λ −ΦM,DȳD ] .
Thus, we have
wM,D,λ − w̄tM,D,λ =w̄t−1M,D,λ − (CM,D + λI) −1 [ (CM,D + λI)w̄ t−1 M,D,λ −ΦM,DȳD ] − w̄t−1M,D,λ +
m∑ j=1 |Dj | |D| (CM,Dj + λI) −1 [ (CM,D + λI)w̄ t−1 M,D,λ −ΦM,DȳD ] =
m∑ j=1 |Dj | |D| [ (CM,Dj + λI) −1 − (CM,D + λI)−1 ] [ (CM,D + λI)w̄ t−1 M,D,λ −ΦM,DȳD ] =
m∑ j=1 |Dj | |D| (CM,Dj + λI) −1 [CM,D − CM,Dj ] (CM,D + λI)−1 [(CM,D + λI)w̄t−1M,D,λ −ΦM,DȳD]
= m∑ j=1 |Dj | |D| (CM,Dj + λI) −1 [CM,D − CM,Dj ] [w̄t−1M,D,λ −wM,D,λ]
= m∑ j=1 |Dj | |D| (CM,Dj + λI) −1 [CM,D − CM ] [w̄t−1M,D,λ −wM,D,λ]
+ m∑ j=1 |Dj | |D| (CM,Dj + λI) −1 [CM − CM,Dj ] [w̄t−1M,D,λ −wM,D,λ]
:= m∑ j=1 |Dj | |D| ℵj1 + ℵ j 2.
(34)
Note that
SMℵj1 = SM (CM + λI)−1/2(CM + λI)1/2(CM,Dj + λI)−1/2
(CM,Dj + λI) −1/2(CM + λI) 1/2(CM + λI) −1/2 [CM,D − CM ] (CM + λI)−1/2
(CM + λI) −1/2(CM + λI) ( w̄t−1M,D,λ −wM,D,λ ) .
(35)
Note that ‖SM (CM + λI)−1/2‖ = ‖(CM + λI)−1/2CM (CM + λI)−1/2‖1/2 ≤ 1, so, we have ‖SMℵj1‖ρ ≤ J 2M,DjQM,D ∥∥∥(CM + λI)−1/2(CM + λI)(w̄t−1M,D,λ −wM,D,λ)∥∥∥ . (36)
Note that CM = S∗MSM , so
CM ( w̄t−1M,D,λ −wM,D,λ ) = S∗MSM ( w̄t−1M,D,λ −wM,D,λ ) = S∗M ( f̄ t−1M,D,λ − fM,D,λ ) .
Substituting the above inequality into Eq. 36, we have ‖SMℵj1‖ρ ≤ J 2M,DjQM,D ∥∥∥(CM + λI)−1/2S∗M (f̄ t−1M,D,λ − fM,D,λ)∥∥∥
+ λJ 2M,DjQM,D ∥∥∥(CM + λI)−1/2 (w̄t−1M,D,λ −wM,D,λ)∥∥∥
≤ J 2M,DjQM,D ( ‖f̄ t−1M,D,λ − fM,D,λ‖ρ + √ λ‖w̄t−1M,D,λ −wM,D,λ‖ ) ,
the last inequality use the fact that ‖(CM + λI)−1/2S∗M‖ = ‖(CM + λI)−1/2CM (CM + λI)−1/2‖1/2 ≤ 1. Using the same process, we can obtain that
‖SMℵj2‖ρ ≤ J 2M,DjQM,Dj ( ‖f̄ t−1M,D,λ − fM,D,λ‖ρ + √ λ‖w̄t−1M,D,λ −wM,D,λ‖ ) .
Thus, we have ‖fM,D,λ − f̄ tM,D,λ‖ρ =‖SM (wM,D,λ − w̄tM,D,λ)‖ρ ≤ m∑ j=1 |Dj | |D| ‖SMℵj1‖ρ + ‖SMℵ j 2‖ρ
≤ m∑ j=1 |Dj | |D|
(J 2M,DjQM,D + J 2 M,DjQM,Dj )(
‖f̄ t−1M,D,λ − fM,D,λ‖ρ + √ λ‖w̄t−1M,D,λ −wM,D,λ‖ ) .
(37)
According to 34, we know that wM,D,λ − w̄tM,D,λ = m∑ j=1 |Dj | |D| (ℵj1 + ℵ j 2)
= m∑ j=1 |Dj | |D| (CM,Dj + λI) −1 [CM,D − CM ] [w̄t−1M,D,λ −wM,D,λ]
+ m∑ j=1 |Dj | |D| (CM,Dj + λI) −1 [CM − CM,Dj ] [w̄t−1M,D,λ −wM,D,λ]
= m∑ j=1 |Dj | |D| (CM,Dj + λI) −1(CM + λI)(CM + λI) −1 [CM,D − CM ] [w̄t−1M,D,λ −wM,D,λ]
+ m∑ j=1 |Dj | |D| (CM,Dj + λI) −1(CM + λI)(CM + λI) −1 [CM − CM,Dj ] [w̄t−1M,D,λ −wM,D,λ]. Thus, one can obtain that
‖wM,D,λ − w̄tM,D,λ‖ ≤ m∑ j=1 J 2M,Dj (LM,D + LM,Dj )‖w̄ t−1 M,D,λ −wM,D,λ‖,
where LM,D := ‖(CM + λI)−1(CM −CM,D)‖. From Proposition 6, we know that LM,D = QM,D. Thus, we have
‖wM,D,λ − w̄tM,D,λ‖ ≤ m∑ j=1 J 2M,Dj (QM,D +QM,Dj )‖w̄ t−1 M,D,λ −wM,D,λ‖. (38)
Combining 37 and 38, we have
‖fM,D,λ − f̄ tM,D,λ‖ρ + √ λ‖wM,D,λ − w̄tM,D,λ‖
≤ m∑ j=1 |Dj | |D| (J 2M,DiQM,D + J 2 M,DjQM,Dj ) ( ‖f̄ t−1M,D,λ − fM,D,λ‖ρ + √ λ‖w̄t−1M,D,λ −wM,D,λ‖ ) +
m∑ j=1 |Dj | |D| J 2M,Dj ( QM,D +QM,Dj )√ λ‖w̄t−1M,D,λ −wM,D,λ‖
≤ m∑ j=1 |Dj | |D| (2J 2M,DjQM,Dj + 2J 2 M,DjQM,Dj ) ( ‖f̄ t−1M,D,λ − fM,D,λ‖ρ + √ λ‖w̄t−1M,D,λ −wM,D,λ‖ )
≤ 2 m∑ j=1 |Dj | |D| J 2M,DjQM,D + J 2 M,DjQM,Dj t (‖f̄0M,D,λ − fM,D,λ‖ρ +√λ‖w̄0M,D,λ −wM,D,λ‖) , which prove the result.

heading: E.2 APPENDIX: PROOF OF THEOREM 6
text: Proof of Theorem 6. Substituting Propositions 4, 5, 7, 8 and Lemma 4 into Proposition 9, we have
‖f̄0M,D,λ − fM,D,λ‖ρ + √ λ‖w̄0D,λ −wM,D,λ‖2
=O  m∑ j=1 |Dj | |D| ( QM,D +QM,Dj ) (RM,D +KM,D + ‖fM,λ − fρ‖ρ)  =O
 m∑ j=1 |Dj | |D| (QM,Dj +QM,Dj ) ( ΥM,Dj ,λ + ‖fM,λ − fλ‖+ ‖fλ − fρ‖ρ ) . Combining the above inequality and Proposition 10, and note that QM,D ≤ QM,Dj , we can obtain that
‖f̄ tM,D,λ − fM,D,λ‖ρ
≤O  m∑ j=1 |Dj | |D| QM,Dj t m∑ j=1 |Dj | |D| QM,Dj ( ΥM,Dj ,λ + ‖fM,λ − fλ‖+ ‖fλ − fρ‖ρ ) . Note that
‖f̄ tM,D,λ − fρ‖ρ = ‖f0M,D,λ − fM,D,λ + fM,D,λ − fM,λ + fM,λ − fλ + fλ − fρ‖ρ ≤‖f̄ tM,D,λ − fM,D,λ‖ρ + ‖fM,D,λ − fM,λ‖ρ + ‖fM,λ − fλ‖ρ + ‖fλ − f‖ρ.
(39)
Thus, by Lemmas 2, 3, one can see that if M ≥ Ω (( N 2r−1(λ) λ2r−1 ) N∞(λ)2−2r ∨N∞(λ) ) ,
we have
‖f̄ tM,D,λ − fρ‖ρ
=O  m∑ j=1 |Dj | |D| QM,Dj t m∑ j=1 |Dj | |D| QM,Dj ( ΥM,Dj ,λ + λ r )+ ΥM,D,λ log 1 δ + λr  . From 33, by setting |D1| = . . . = |Dm|, λ = Ω(|D|− 1 2r+γ ), we know that
ΥM,Dj ,λ = O (√ m|D|− r 2r+γ +m|D|− 4r+2γ−1 4r+2γ ) ,
QM,Dj = O ( m|D|− 2r+γ−1 2r+γ + √ m|D|− 2r+γ−1 4r+2γ ) ,
ΥM,D,λ = O ( |D|− r 2r+γ + |D|− 4r+2γ−1 4r+2γ ) .
Thus, when m ≤ Ω ( |D| (2r+γ−1)(t+1) (2r+γ)(t+2) ) , we have ‖f̄ tM,D,λ − fρ‖ρ = O ( |D|− r 2r+γ ) , which proves the result.

heading: F APPENDIX: PROOF OF THEOREMS 1, 2, 3
text: Proof. From (Smale & Zhou, 2007; Caponnetto & Vito, 2007), if r = 1/2, then fρ ∈ HK . Thus, in this case, fHK exists and E(fHK ) = E(fρ). Note that Assumption 1 is always satisfied for γ = 1 and c = τ2. So, using Theorems 4, 5 and 6 with r = 1/2, γ = 1 and c = τ2, Theorem 1, 2, and 3 can be proved.

Review summaries:
The paper analyses generalization properties of distributed kernel ridge regression (DKRR) with random features and communications. It studies optimal learning rates of the generalization bounds both in expectation and in probability. In the case of DKRR with random features, the optimal learning rate in expectation is shown to achieve by relaxing the requirement on the number of partitions from $O(1)$ (Li et al., 2019a) to $O(|D|^{0.5})$ (Theorem 1). Within the same setup of random features, the number of partitions is relaxed to $O(|D|^{0.25})$ guaranteeing optimal generalization performance in probability (Theorem 2). The latter bound $O(|D|^{0.25})$ on partition count is much smaller then $O(|D|^{0.5})$. However, as proved in Theorem 3, allowing multiple communication rounds in DKRR-RF, up to $O(|D|^{0.5})$ partitions can be handled depending on the number of communication rounds. In other words, it can exploit more partitions at the cost of more communication rounds.

The paper investigates an algorithm for distributed learning with random Fourier features. The main idea is to sample M random Fourier features and split the data into m chunks. Each chunk is processed on a separate machine that outputs a linear hypothesis using the sampled M random features. The hypotheses coming from different machines are then aggregated on the master machine via importance weighting. In particular, each hypothesis is assigned importance weight proportional to its data chunk size (see Eq. 3). The regularization parameter is fixed across different machines. The main contribution of the work is a consistency bound. In comparison to a previous bound on the divide & conquer algorithm (Li et al., arXiv 2019), this one does not require a constant number of machines (in my understanding of the related work section).

This paper studies the statistical properties of distributed kernel ridge regression together with random features (DKRR-RF), and obtain optimal generalization bounds under the basic setting in the attainable cases.  Numerical results are given for the studied new algorithms. The algorithms and the derived results are new and interesting to me. However, the presentations as well as the citations need some major revision before the publication.