Paper deep dive
Improving TensorSketch Using Complex Random Variables
Amit Sharma, Mohammad Azhar Khan, Rameshwar Pratap, Keegan Kang
Intelligence
Status: succeeded | Model: Gemma-4-26B-A4B | Prompt: intel-v1 | Confidence: 93%
Last extracted: 8/16/2026, 2:43:22 AM
Summary
This paper introduces a Complex-to-Real (CtR) variant of the TensorSketch algorithm for approximating high-dimensional polynomial kernels. While previous complex-valued sketching methods improved variance bounds from $3^p/D$ to $2^p/D$, they relied on dense Johnson-Lindenstrauss projections with high computational cost. The proposed method integrates complex random variables (fourth roots of unity) into the sparse TensorSketch framework, achieving the improved $2^p/D$ variance bound while retaining the input-sparsity running time $O(p( ext{nnz}(oldsymbol{x}) + D ext{log} D))$. Theoretical proofs and experiments on synthetic and real-world datasets validate the accuracy and efficiency of the new estimator.
Entities (8)
Relation Signals (8)
Complex-to-Real TensorSketch → isvariantof → TensorSketch
confidence 98% · In this work, we introduce a simple variant of TensorSketch
Complex-to-Real TensorSketch → improves → Variance Bound
confidence 95% · The proposed construction obtain variance improvements comparable to the complex JL-type estimator... reducing it from 3^p to 2^p.
Complex-to-Real TensorSketch → retains → Input-Sparsity Running Time
confidence 95% · while retaining its advantage of the input-sparsity running time.
Complex-to-Real TensorSketch → targets → Polynomial Kernel
confidence 95% · TensorSketch by pham2013fast, kar2012random provides efficient sketching algorithms for high-dimensional polynomial kernels
Complex-to-Real TensorSketch → utilizes → Complex Random Variables
confidence 94% · It combines a random function whose values are drawn independently and uniformly from the fourth roots of unity
TensorSketch → uses → CountSketch
confidence 92% · pham2013fast extends the sparse CountSketch algorithm
Johnson-Lindenstrauss Transform → hascomputationalcost → O(pDd)
confidence 90% · kar2012random uses dense Johnson-Lindenstrauss (JL)-type projections with computational cost O(pDd)
TensorSketch → hascomputationalcost →
Cypher Suggestions (0)
No Cypher suggestions yet.
Abstract
Abstract:\texttt{TensorSketch} by~\cite{pham2013fast,kar2012random} provides efficient sketching algorithms for high-dimensional polynomial kernels $\vec{x}^{\otimes p} \in \R^{d^p}$. \cite{kar2012random} uses dense Johnson-Lindenstrauss (JL)-type projections with computational cost $O(pDd)$, where $D$ denotes the sketch dimension, whereas~\cite{pham2013fast} extends the sparse \texttt{CountSketch}~\citep{count_sketch} algorithm, yielding a faster algorithm for high-dimensional sparse inputs with running time $O\big(p(\nnz{\vec{x}} + D \log D)\big)$. However, the variance of both estimators grows exponentially with the polynomial degree $p$, scaling as $3^{p}/D$. Recent work by~\cite{pmlr-v206-wacker23a} showed that using complex-valued distribution reduces this dependence to $2^{p}/D$ for the approach of~\cite{kar2012random}. However, their method relies on dense JL-type projections with computational cost $O(pDd)$ and does not extend to the algorithm of~\cite{pham2013fast}. In this work, we introduce a simple variant of \texttt{TensorSketch}~\citep{pham2013fast} that achieves the same variance bound as~\cite{pmlr-v206-wacker23a}, while retaining its advantage of the input-sparsity running time. We validate our results with supporting experiments on synthetic and real-world datasets.
Tags
Links
- Source: https://arxiv.org/abs/2608.10523v1
- Canonical: https://arxiv.org/abs/2608.10523v1
Trouble viewing inline? Open PDF directly →
Full Text
97,757 characters extracted from source content.
Expand or collapse full text
Improving TensorSketch Using Complex Random Variables Amit Sharma* Department of Computer Science and Engineering IIT Hyderabad India Mohammad Azhar Khan* Department of Computer Science and Engineering IIT Hyderabad India Rameshwar Pratap Department of Computer Science and Engineering IIT Hyderabad India Keegan Kang Department of Mathematics and Statistics Bucknell University USA Abstract TensorSketch by pham2013fast, kar2012random provides efficient sketching algorithms for high-dimensional polynomial kernels ⊗p∈ℝdpx p ^d^p. kar2012random uses dense Johnson-Lindenstrauss (JL)-type projections with computational cost O(pDd)O(pDd) , where D denotes the sketch dimension, whereas pham2013fast extends the sparse CountSketch [count_sketch] algorithm, yielding a faster algorithm for high-dimensional sparse inputs with running time O(p(nnz()+DlogD))O (p(nnz (x )+D D) ). However, the variance of both estimators grows exponentially with the polynomial degree p, scaling as 3p/D3^p/D. Recent work by pmlr-v206-wacker23a showed that using complex-valued distribution reduces this dependence to 2p/D2^p/D for the approach of kar2012random. However, their method relies on dense JL-type projections with computational cost O(pDd)O(pDd) and does not extend to the algorithm of pham2013fast. In this work, we introduce a simple variant of TensorSketch [pham2013fast] that achieves the same variance bound as pmlr-v206-wacker23a, while retaining its advantage of the input-sparsity running time. We validate our results with supporting experiments on synthetic and real-world datasets. †footnotetext: * Equal contribution. 1 Introduction Polynomial kernels (Definition 3) are widely used in machine learning to model higher-order, non-linear interactions among input features. For vectors ,∈ℝdx,y ^d, a degree-p polynomial kernel is defined as k(,)=⟨,⟩pk(x,y)= ,y ^p. This kernel is equivalent to mapping ∈ℝdx ^d to its p-fold Kronecker product ⊗p∈ℝdp,x p ^d^p, which is ⟨⊗p,⊗p⟩≡⟨,⟩p. p,y p ≡ ,y ^p. As the dimension dpd^p grows exponentially with p, directly computing these inner products is computationally infeasible. A dense Johnson-Lindenstrauss (JL) transform [johnson1984extensions] can reduce their dimensionality while approximately preserving their pairwise inner product, however applying it requires time proportional to dpd^p, which is exponential in p. To address this inefficiency, prior works such as kar2012random and pham2013fast proposed randomized sketching techniques to approximate polynomial kernels efficiently. kar2012random introduced random feature maps based on JL-type projections. They define a randomized linear map :ℝdp→ℝDS:R^d^p ^D by (⊗p):=(1⊙⋯⊙p)/DS(x p):= (W_1x\, \,·s\, \,W_px )/ D, where each i∈ℝD×dW_i ^D× d is a random projection matrix (e.g., Gaussian or Rademacher with i.i.d. entries), and ⊙ denotes the element-wise (Hadamard) product. Such that k^(,):=⟨(⊗p),(⊗p)⟩. k (x,y ):= (x p),S(y p) . This estimator can be termed as a JL-type variant of TensorSketch, and its computational cost scales with O(pDd)O(pDd). pham2013fast propose TensorSketch, which combines CountSketch [count_sketch] with Fast Fourier Transform (FFT)-based convolution to compute the sketch implicitly. This avoids forming the full vector and runs in input-sparsity time O(p(nnz()+DlogD)),O (p(nnz (x )+D D) ), making it preferred over kar2012random for high-dimensional sparse data. The variance of both algorithms kar2012random, pham2013fast, pham2025tensorsketchfastscalable grows as 3p/D3^p/D. Recent progress by wacker2022improved introduced a complex-valued variant of JL-type TensorSketch. It improves the variance dependence of the sketch on the polynomial degree, reducing it from 3p3^p to 2p2^p. This improvement is achieved by incorporating complex-valued randomness into the sketch. However, the resulting sketch vectors are complex-valued and cannot be directly compared with their real-valued counterparts. To address this, pmlr-v206-wacker23a proposed the Complex-to-Real (CtR) construction. It first computes a complex random feature map of the embedding dimension half of the size of its real counterpart and then forms a real embedding by concatenating its real and imaginary parts. This preserves inner products and achieves improved variance bounds while yielding a real-valued sketch. However, as a JL-type method, it still requires dense multiplications with cost O(pDd)O(pDd), which can be inefficient for high-dimensional sparse data. In this work, we address the above limitation by developing the Complex-to-Real (CtR) variant of TensorSketch. The proposed construction obtain variance improvements comparable to the complex JL-type estimator of pmlr-v206-wacker23a, while retaining the input-sparsity running time guarantees of TensorSketch [pham2013fast]. Contributions. We propose a Complex-to-Real variant of TensorSketch for degree-p polynomial kernels. It combines a random function whose values are drawn independently and uniformly from the fourth roots of unity with FFT-based tensor sketching while producing real-valued embeddings. We prove that the resulting estimator is unbiased for ⟨⊗p,⊗p⟩ p,y p . We also derive an upper bound on variance with 2p2^p-type dependence on the degree, improving over the classical TensorSketch [pham2013fast] variance scaling, while retaining the sketching time O(p(nnz()+DlogD))O(p(nnz (x )+D D)). The algorithm is defined in Definition 8, and its theoretical guarantees are stated in Theorem 2. We note that the use of complex-valued random variables in sketching algorithms has been explored in prior work. For example, pmlr-v206-wacker23a, wacker2022improved employ complex random variables to construct sketches for polynomial kernels, while meyer2026hutchinson use them for trace estimation of implicit matrices. The key idea in pmlr-v206-wacker23a, wacker2022improved is to sketch ⊗px p using element-wise products of independent sketches of the factors (1,…,p)(x_1,…,x_p). Owing to this independence structure, Khintchine’s inequality (Definition 4) can be applied, and a simple induction on p yields a tighter upper bound in the complex setting than in the real-valued case. In contrast, TensorSketch is built upon the hash-based CountSketch framework, where the sketch components are not independent and no direct analogue of Khintchine’s inequality is available. As a result, the techniques from prior work do not extend to our setting. Moreover, our analysis of CountSketch with complex-valued random variables (Theorem 5, Appendix) shows that, by itself, the complex construction of CountSketch provides no variance reduction over its real-valued version. Therefore, our result demonstrating improved variance for TensorSketch with complex random variables is non-trivial and stems from the vanishing of certain cross terms in the variance analysis, which does not occur in the real-valued setting. Implications of our results. TensorSketch [pham2013fast] enables efficient training of linear SVMs [sun2018but, li2019towards], supports scalable deep learning and the theoretical analysis of over-parameterized neural networks [yehudai2019power, zandieh2021scaling], and has been successfully applied to compact bilinear pooling for fine-grained visual recognition [gao2016compact] as well as multimodal fusion architectures [fukui2016multimodal]. Polynomial kernels themselves are widely adopted in applications including natural language processing [goldberg2008splitting], recommender systems [rendle2010factorization], and genomics [aschard2016incorporating]. Our proposed CtR TensorSketch achieves more accurate estimates than TensorSketch [pham2013fast] while maintaining the same asymptotic time complexity, making it a promising alternative for the above applications. Organization of the paper. The remainder of this paper is organized as follows. Section 2 reviews prior work on polynomial kernel approximation and randomized sketching. Section 3 presents necessary background, including CountSketch, TensorSketch, polynomial kernels, and the CtR framework. Section 4, proposes our algorithm and states theoretical guarantees on its accuracy and efficiency. Section 5 reports empirical comparisons of real and JL-type methods with our approach. Section 6 concludes and outlines some future research directions. 2 Related Work Polynomial kernels can be expressed via tensor feature maps, where the representation is the tensor product ⨂i=1pi _i=1^px_i of vectors 1∈ℝd1,…,p∈ℝdpx_1 ^d_1,…,x_p ^d_p to capture higher–order interactions. Explicit construction requires O(∏i=1pdi)O\! ( _i=1^pd_i ) memory and is infeasible even for moderate p or dimensions d. To avoid this blow-up, prior work proposes sketching methods that compute (⨂i=1pi)S\! ( _i=1^px_i ) without forming the full tensor, including JL-type product embeddings [kar2012random] and hashing-based TensorSketch constructions built on CountSketch [pham2013fast, pham2025tensorsketchfastscalable]. For ,∈ℝdx,y ^d, both TensorSketch estimators have variance ≤3p−1D(‖22‖22),≤ 3^p-1D\! ( \|x \|_2^2 \|y \|_2^2 ), which yields 3p3^p-type variance growth. Both JL-type and hashing-based estimators scale as Θ(3p/D) (3^p/D) for degree-p features, but TensorSketch [pham2013fast, pham2025tensorsketchfastscalable] is more efficient due its input-sparsity time complexity. Recent work shows that complex-valued distributions reduce the variance of randomized feature maps. Complex JL-type constructions for polynomial kernels achieve lower variance than their real-valued counterparts [wacker2022improved]. In particular, replacing real Rademacher or Gaussian variables with complex-valued distributions improves the variance dependence from 3p3^p-type to 2p2^p-type. Further work by pmlr-v206-wacker23a introduced the Complex-to-Real (CtR) construction, preserving the variance improvement of complex-valued sketches while yielding real-valued embeddings. These embeddings enables direct comparison with real sketches. However, these approaches rely on dense JL-style projections and therefore incur higher computational cost especially for high-dimensional sparse data. Building on the hashing-based sketching framework of [pham2013fast], we introduce a Complex-to-Real (CtR) variant of TensorSketch for polynomial kernel approximation. Our method uses fourth roots-of-unity instead of real Rademacher variables and applies a structured complex-to-real transformation, yielding unbiased real-valued embeddings with improved variance dependence. In contrast, our method achieves variance reduction through a lightweight modification of the sketch construction. This is followed by a Complex-to-Real (CtR) conversion to produce a real-valued embedding that makes it possible to compare it fairly with its real counterpart. The resulting sketch preserves the original structure and input-sparsity running time while achieving provable variance improvements over original TensorSketch [pham2013fast]. Variance reduction for randomized sketching has been widely studied using statistical variance reduction techniques such as control variates (CV) and maximum likelihood estimation (MLE). The CV method reduces variance by leveraging a correlated auxiliary variable with a known expectation, while MLE estimates unknown parameters by maximizing the likelihood of the observed data, often yielding statistically efficient estimators. These techniques have been extensively investigated for a variety of classical randomized sketching methods, including Johnson–Lindenstrauss (JL) transforms Li2006very, CountSketch DBLP:conf/uai/PratapK21, Tug-of-War sketches DBLP:conf/acda/PratapVK21, signed random projections Pmlr-v80-kang18b, feature hashing DBLP:journals/ml/VermaPT22, and compressed matrix multiplication DBLP:journals/ipl/VermaDPT25, among others. However, to the best of our knowledge, analogous variance reduction techniques have not yet been developed for TensorSketch. While these approaches can substantially reduce estimator variance, they typically incur additional computational and algorithmic overhead. 3 Preliminaries Notation: We use the following notation in the paper. Bold lowercase letters (e.g., ,x,y) denote vectors, bold uppercase letters (e.g., ,S,T) denote matrices. For a positive integer D, we write [D]:=1,2,…,D[D]:=\1,2,…,D\ for the corresponding index set. For d,p∈ℕd,p , we denote by [d]p[d]^p the set of all p-tuples (i1,…,ip)(i_1,…,i_p) with ij∈[d]i_j∈[d]. The sets ℝR, ℤZ and ℂC denote the real, integer and complex domains, respectively. For any vector x, nnz()nnz (x ) denotes the number of its nonzero entries. The symbol ⟨,⟩ ,y denotes the standard inner product, ‖2 \|x \|_2 the Euclidean norm, and ‖F \|S \|_F the Frobenius norm of a matrix. For a complex number z∈ℂz written as z=a+ibz=a+ib, its complex conjugate is defined as z¯ z. The norm of z is |z|=a2+b2=zz¯.|z|= a^2+b^2= z\, z. The Kronecker product is written as ⊗ , while “⋅·” denotes standard matrix multiplication and ⊙ denotes the elementwise (Hadamard) product. The indicator function is denoted by [⋅]1 [· ]. The probabilistic quantities use Pr[⋅]Pr[·] for probability, [⋅]E [· ] for expectation and Var(⋅)Var (· ) for variance. Definition 1 (CountSketch [count_sketch]). Given an input vector ∈ℝdy ^d, the CountSketch is a randomized linear map ∈ℝD×dT ^D× d that maps y to a lower-dimensional vector =∈ℝD.z=Ty ^D. The CountSketch matrix T is constructed by two hash functions: (a) h:[d]→[D]h [d]→[D] a 22-wise independent hash function, and (b) s:[d]→1,−1s:[d]→\1,-1\ a 44-wise independent random sign function. The jthj^th entry of vector ∈ℝDz ^D is computed as, zj=∑h(i)=js(i)yi,∀j∈[D].z_j= _h(i)=js(i)\,y_i,\;∀ j∈[D]. The time complexity of computing the CountSketch is O(nnz())O\! (nnz (y ) ). For pairwise inner–product analysis, let ,∈ℝdx,y ^d and define the CountSketch inner–product estimator as k^(,):=⟨,⟩. k (x,y ):= ,Ty . Then the estimator is unbiased [k^(,)]=⟨,⟩,E [ k (x,y ) ]= ,y , and its variance satisfies Var(k^(,))≤1D(⟨,⟩2+‖22‖22−2∑ixi2yi2).Var ( k (x,y ) )≤ 1D ( ,y ^2+ \|x \|_2^2 \|y \|_2^2-2 _ix_i^2y_i^2 ). Definition 2 (TensorSketch of Degree p [pham2013fast, pham2025tensorsketchfastscalable]). Let p≥2p≥ 2. For each t∈[p]t∈[p], let ht:[d]→[D]h_t:[d]→[D] be 22-wise independent hash functions and σt:[d]→−1,+1 _t:[d]→\-1,+1\ be 44-wise independent random sign functions. The degree-p TensorSketch matrix ∈ℝD×dpS ^D× d^p is defined for all r∈[D]r∈[D] and (i1,…,ip)∈[d]p(i_1,…,i_p)∈[d]^p by r,(i1,…,ip)=(∏t=1pσt(it))[∑t=1pht(it)≡r(modD)].S_r,(i_1,…,i_p)= ( _t=1^p _t(i_t) )1 [ _t=1^ph_t(i_t)≡ r D ]. For any ∈ℝdx ^d, the sketch (⊗p)S\! (x p ) can be computed implicitly in time O(p(nnz()+DlogD))O(p(nnz (x )+D D)) using FFT-based convolution, without explicitly forming ⊗px p. Let ,∈ℝdx,y ^d, then estimate for pairwise inner–product is defined as follows: k^(,):=⟨(⊗p),(⊗p)⟩. k (x,y ):= \! (x p ),S\! (y p ) . Then TensorSketch provides an unbiased estimator of the degree–p polynomial kernel [k^(,)]=⟨,⟩p,E [ k (x,y ) ]= ,y ^p, and its variance satisfies Var(k^(,)) ( k (x,y ) ) ≤1D((2⟨,⟩2+∥22∥22⋯ ≤ 1D ( (2 ,y ^2+ \|x \|_2^2\, \|y \|_2^2·s ⋯−2∑i=1dxi2yi2)p−⟨,⟩2p), ·s-2 _i=1^dx_i^2y_i^2 )^p- ,y ^2p ), ≤3p−1D‖22p‖22p. ≤ 3^p-1D\, \|x \|_2^2p\, \|y \|_2^2p. Definition 3 (Polynomial Kernel [10.7551/mitpress/4175.001.0001]). We consider polynomial kernels of degree p∈ℕp of the form k(,)=(γ⊤+ν)p,k(x,y)= (γ\,x y+ν )^p, for ,∈ℝdx,y ^d with γ,ν≥0γ,ν≥ 0. The parameters γ and ν can be absorbed into the inputs by introducing the augmented vectors ~=(γ⊤,ν)⊤,~=(γ⊤,ν)⊤∈ℝd+1. x= ( γ\,x ,\; ν ) , y= ( γ\,y ,\; ν ) ^d+1. With this transformation, the kernel admits a homogeneous representation, (γ⊤+ν)p=(~⊤~)p=(~⊗p)⊤~⊗p, (γ\,x y+ν )^p=( x y)^p= ( x p ) y p, where ~⊗p x p denotes the p-fold tensor product of ~ x. Hence, without loss of generality, we may treat the polynomial kernel as a homogeneous kernel in the augmented space ℝd+1R^d+1. Definition 4 (Khintchine Inequality haagerup2007best, pmlr-v206-wacker23a). Let =(x1,…,xd)∈ℝdx=(x_1,…,x_d) ^d and let εii=1d\ _i\_i=1^d be independent Rademacher random variables. For every p>0p>0, there exists a constant KpK_p such that (|∑i=1dxiεi|p)1/p≤Kp‖2. (E | _i=1^dx_i _i |^p )^1/p≤ K_p\|x\|_2. (1) The optimal constants are known explicitly. • If εi∈−1,+1 _i∈\-1,+1\ are real-valued Rademacher variables, then Kp=1,0<p≤2,2π−12pΓ(p+12)1p,p>2. K_p= cases1,&0<p≤ 2,\\[5.69054pt] 2\,π^- 12p \! ( p+12 ) 1p,&p>2. cases (2) • If εi _i are complex Rademacher variables taking values in 1,−1,i,−i\1,-1,i,-i\ uniformly at random, then Kp^=Γ(p2+1)1p. K_p= \! ( p2+1 ) 1p. (3) Moreover, for every p>2p>2, the complex constant Kp K_p is strictly smaller than the corresponding real constant KpK_p. Definition 5 (Framework for Complex-to-Real (CtR) Sketches). pmlr-v206-wacker23a suggests a framework of analysing sketching algorithms involving complex random variables. Let z=a+ibz=a+ib be a complex random variable with a,b∈ℝa,b , we have |z|2=a2+b2andRez2=a2−b2. |z |^2=a^2+b^2 \z^2 \=a^2-b^2. Combining both gives a2=12(|z|2+Rez2).a^2= 12 ( |z |^2+Re \z^2 \ ). The scalar a is real-valued and its variance Var(a)=[a2]−[a]2Var (a )=E [a^2 ]-E [a ]^2 is therefore Var(a) (a ) =12Re[|z|2]+[z2]−2[a]2. = 12Re\E[|z|^2]+E[z^2]-2E[a]^2\. (4) Let ,∈ℝdx,y ^d and let ΦC:ℝdp→ℂD2 _C:R^d^p D2. Then k^C(,):=ΦC(⊗p)ΦC(⊗p)¯⊤∈ℂ k_C (x,y ):= _C(x p) _C(y p) is a complex-valued estimate of the polynomial kernel k(,)k(x,y), and hence k^C(,)=Rek^C(,)+iImk^C(,) k_C (x,y )=Re\ k_C (x,y )\+iIm\ k_C (x,y )\. The Complex-to-Real sketch is defined as follows: k^CtR(,):=Rek^C(,), k_CtR (x,y ):=Re\ k_C (x,y )\, =ReΦC(⊗p)⊤ReΦC(⊗p)+ =\!Re \ _C(x p) \ \!Re \ _C(y p) \\!+ ⋯+ImΦC(⊗p)⊤ImΦC(⊗p), ·s+\!Im \ _C(x p) \ \!Im \ _C(y p) \, =ΦCtR()⊤ΦCtR(). = _CtR(x) _CtR(y). where, ΦCtR():=[ReΦC,ImΦC]∈ℝD. _CtR(x):=[Re \ _C \,Im \ _C \] ^D. We now derive the variance of k^CtR(,) k_CtR (x,y ) using Equation (4): Var(k^CtR(,)) ( k_CtR (x,y ) ) =12Re[|k^C|2]+[k^C2]+2[Rek^C]2 = 12Re\E[| k_C|^2]+E[ k_C^2]+2E[Re\ k_C\]^2\ =12Re[|k^C|2]+[k^C2]+2[k^C]2, = 12Re\E[| k_C|^2]+E[ k_C^2]+2E[ k_C]^2\, (5) where, k^C:=k^C(,)∈ℂ. k_C\!:=\! k_C (x,y )\!∈\!C. This framework plays a key role in computing the expectation and variance of k^CtR(,) k_CtR (x,y ). Definition 6 (Complex-to-Real Polynomial Sketch [pmlr-v206-wacker23a]). Let p∈ℕp and D=2kD=2k for some k∈ℕk . For each i∈[p]i∈[p], let i∈ℂD/2×dW_i ^D/2× d be a random matrix whose rows are sampled independently from a zero-mean distribution satisfying [∗]=dE [ww^* ]=I_d (e.g., complex Gaussian or complex Rademacher) where ∗* refers to conjugate-transpose. Define the complex random feature map ΦC(⊗p):=2D(⊙⊙⋯⊙)∈ℂD/2, _C(x p):= 2D\; (W_1x _2x ·s _px )\;∈\;C^D/2, Then, the Complex-to-Real (CtR) sketch of ∈ℝdpx ^d^p is the real-valued vector ΦCtR(⊗p)∈ℝD _CtR(x p) ^D defined by ΦCtR(⊗p):=[ReΦCtR(⊗p)ImΦCtR(⊗p)]∈ℝD. _CtR(x p):= bmatrixRe \ _CtR(x p) \\\ Im \ _CtR(x p) \ bmatrix ^D. For any ⊗p,⊗p∈ℝdpx p,y p ^d^p, the kernel estimator can be defined as k^CtR(,) k_CtR (x,y ) =ΦCtR()⊤ΦCtR(). = _CtR(x) _CtR(y). For both Gaussian and Rademacher constructions, the CtR estimator satisfied the following: [k^CtR(,)] [ k_CtR (x,y ) ] =⟨,⟩p, = ,y ^p, Var(k^CtR(,)) ( k_CtR (x,y ) ) ≤2p+1−2D‖22p‖22p. ≤ 2^p+1-2D\,\|x\|_2^2p\,\|y\|_2^2p. 4 Complex-to-Real(CtR) TensorSketch Algorithm Variance Sketching Time TensorSketch(CtR) 2p+1−2D‖22p‖22p 2^p+1-2D \|x \|_2^2p\, \|y \|_2^2p O(p(nnz()+nnz()+DlogD))O(p(nnz (x )+nnz (y )+D D)) TensorSketch(Real)[pham2013fast, pham2025tensorsketchfastscalable] 3p−1D‖22p‖22p 3^p-1D \|x \|_2^2p\, \|y \|_2^2p O(p(nnz()+nnz()+DlogD))O(p(nnz (x )+nnz (y )+D D)) JL(CtR Radamacher)[pmlr-v206-wacker23a] 2p+1−2D‖22p‖22p 2^p+1-2D \|x \|_2^2p\, \|y \|_2^2p O(pDd)O(pDd) JL(CtR Guassian)[pmlr-v206-wacker23a] 2p+1−2D‖22p‖22p 2^p+1-2D \|x \|_2^2p\, \|y \|_2^2p O(pDd)O(pDd) Table 1: Comparison of variance bounds and sketching time for CtR TensorSketch and baseline methods. Here, ,∈ℝdx,y ^d are input vectors, p denotes the polynomial degree, and D is the sketch dimension. All methods provide unbiased estimates of the inner product between the vectors ⊗p,⊗p∈ℝdpx p,y p ^d^p, i.e., ⟨⊗p,⊗p⟩ p,y p . In this section, we first define Complex CountSketch (Definition 7). Building on its circular convolution structure, we then introduce Complex-to-Real TensorSketch. (Definition 8) along with the associated kernel. Theorem 2 discusses the guarantees of the estimators proposed in Definition 8. Definition 7 (Complex CountSketch Mapping). Let ∈ℝdx ^d and D be the sketching dimension. Sample a Complex CountSketch matrix ∈ℂD×dC ^D× d with the following two independent functions • h:[d]→[D]h:[d]→[D] is an universal hash function that assigns each coordinate independently and uniformly to one of the D buckets, and • s:[d]→1,ω,ω2,ω3s:[d]→\1,ω,ω^2,ω^3\ is a random function whose values are drawn independently and uniformly from the four fourth roots of unity. The j-th coordinate of the sketched vector Cx is given by ()j=∑i=1ds(i) 1[h(i)=j]xi,∀j∈[D]. (Cx)_j= _i=1^ds(i)\,1 [h(i)=j ]x_i\,, ∀\,j∈[D]. We define CtR TensorSketch, which extends the Complex CountSketch by applying p independent sketches and combining them using FFT-based convolution to efficiently sketch the vector ⊗p,⊗p∈ℝdpx p,y p ^d^p. Definition 8. (Complex-to-Real (CtR) TensorSketch Mapping) Let ∈ℝdx ^d and fix an integer p≥1p≥ 1. Define a Complex TensorSketch mapping :ℝdp→ℂD/2C:R^d^p ^D/2 constructed from p independent Complex CountSketch mappings. For each r∈[p]r∈[p], let • hr:[d]→[D/2]h_r:[d]→[D/2] be an universal hash function that assigns each coordinate independently and uniformly to one of the D buckets, and • sr:[d]→1,ω,ω2,ω3s_r:[d]→\1,ω,ω^2,ω^3\ be a random function whose values are drawn independently and uniformly from the fourth roots of unity. For each r∈[p]r∈[p], let r∈ℂD/2×dC_r ^D/2× d denote the Complex CountSketch matrix induced by functions (hr,sr)(h_r,s_r) as in Definition 7. The TensorSketch of ⊗p∈ℝdpx p ^d^p is defined implicitly via convolution of the p sketches, and can be computed efficiently using the Fast Fourier Transform as ΦC(⊗p) _C(x p) :=⊗p, :=Cx p, (6) =FFT−1(⨀r=1pFFT(r))∈ℂD/2, =FFT^-1\! ( _r=1^pFFT\! (C_rx ) ) ^D/2, (7) where the product is taken element-wise. We define the CtR sketch as ΦCtR(⊗p):=(ReΦC(⊗p)1,…,ReΦC(⊗p)D/2, \! _CtR(x p)\!:=\! (Re\ _C(x p)_1\,…,Re\ _C(x p)_D/2\, ImΦC(⊗p)1,…,ImΦC(⊗p)D/2)⊤∈ℝD. \ _C(x p)_1\,…,Im\ _C(x p)_D/2\ ) ^D. A detailed proof of results (Lemma 1, 4, and Theorem 2) stated in this section is presented in Appendix A.1. Lemma 1. Let ω=e2πi/4=iω=e^2π i/4=i, and let s:[d]→1,ω,ω2,ω3s:[d]→\1,\ ω,\ ω^2,\ ω^3\ be a random function such that the values s(i)i∈[d]\s(i)\_i∈[d] are drawn independently and uniformly from the four fourth roots of unity. Then, for every i∈[d]i∈[d], the following identities hold [s(i)] [s(i) ] =0, =0, [|s(i)|2] [ |s(i) |^2] =1, =1, [s(i)2] [s(i)^2 ] =0. =0. Furthermore, for any pair of distinct indices i≠ji≠ j, independence implies [s(i)s(j)¯]=[s(i)][s(j)¯]=0. [s(i)\, s(j)]=E[s(i)]\,E [ s(j) ]=0. Theorem 2 establishes the unbiasedness, variance bound, and sketching time of the CtR TensorSketch. Theorem 2. Let ,∈ℝdx,y ^d and ⊗p,⊗p∈ℝdpx p,y p ^d^p. Let ΦCtR _CtR denote a Complex-to-Real TensorSketch as stated in Definition 8. Let k^CtR(,):=ΦCtR(⊗p)⊤ΦCtR(⊗p) k_CtR (x,y ):= _CtR(x p) _CtR(y p), then k^CtR(,) k_CtR (x,y ) satisfies [k^CtR(,)] [ k_CtR (x,y ) ] =⟨⊗p,⊗p⟩=⟨,⟩p, = p,y p = ,y ^p, Var(k^CtR(,)) ( k_CtR (x,y ) ) ≤2p+1−2D‖22p‖22p. ≤ 2^p+1-2D\, \|x \|_2^2p\, \|y \|_2^2p. Moreover, the sketches ΦCtR(⊗p) _CtR(x p) and ΦCtR(⊗p) _CtR(y p) can be computed in O(p(nnz()+DlogD))andO(p(nnz()+DlogD))O(p(nnz (x )+D D))\ and\,O(p(nnz (y )+D D)) time, respectively. Proof. According to Definition 5, we have k^CtR(,):=Rek^C(,) k_CtR (x,y ):=Re \ k_C (x,y ) \, and therefore we first analyze k^C(,) k_C (x,y ). Recall that k^C(,):=ΦC(⊗p)ΦC(⊗p)¯⊤ k_C (x,y ):= _C(x p) _C(y p) , where ΦC(⊗p):=⊗p _C(x p):=Cx p and ΦC(⊗p):=⊗p _C(y p):=Cy p according to Definition 8. In particular, ⊗p,⊗p∈ℂD/2Cx p,Cy p ^D/2 denote the complex CountSketch of the vectors ⊗p,⊗p∈ℝdpx p,y p ^d^p, constructed using the derived functions H:[d]p↦[D/2]H:[d]^p [D/2] and S:[d]p→1,ω,ω2,ω3S:[d]^p→\1,ω,ω^2,ω^3\ defined as H(i1,…,ip) H(i_1,…,i_p) =(∑j=1phj(ij))modD/2, = ( _j=1^ph_j(i_j) ) D/2, S(i1,…,ip) S(i_1,…,i_p) =∏j=1psj(ij). = _j=1^ps_j(i_j). In the proof, we denote X:=⊗pX:=x p and Y:=⊗pY:=y p. Let u,v∈[d]pu,v∈[d]^p be the multi-indices corresponding to the entries of X,Y∈ℝdpX,Y ^d^p. For a multi-index u=(j1,…,jp)u=(j_1,…,j_p), the entry of X is defined as Xu=∏k=1pxjkX_u= _k=1^px_j_k, where the associated linearized index is u=1+∑k=1p(jk−1)∏k′=1p−kdu=1+ _k=1^p(j_k-1) _k =1^p-kd. The entries of Y are defined analogously. We begin by analyzing k^C(,) k_C (x,y ) as defined below, k^C(,):=ΦC(X)ΦC(Y)¯⊤=⟨X,Y¯⟩, k_C (x,y ):= _C(X) _C(Y) = , CY , =∑u,v∈[d]pXuYvS(u)S(v)¯ 1[H(u)=H(v)], = _u,v∈[d]^pX_u\,Y_v\,S(u)\, S(v)\ 1 [H(u)=H(v) ], =⟨X,Y⟩+∑u≠vXuYvS(u)S(v)¯ 1[H(u)=H(v)]. =\ X,Y \!+\! _u≠ vX_u\,Y_v\,S(u)\, S(v)\,1 [H(u)=H(v) ]. Using Lemma 1, we have [S(u)S(v)¯]=0E [S(u)\, S(v) ]=0 for all u≠vu≠ v. Hence, [k^C(,)]=⟨X,Y⟩=⟨⊗p,⊗p⟩=⟨,⟩p. [ k_C (x,y ) ]= X,Y = p,y p = ,y ^p. According to Definition 5, we know k^CtR(,):=Rek^C(,) k_CtR (x,y ):=Re \ k_C (x,y ) \, we have [k^CtR(,)]=[Rek^C(,)]=⟨,⟩p. [ k_CtR (x,y ) ]=E [Re \ k_C (x,y ) \ ]= ,y ^p. (8) Using Equation 5, the variance of k^CtR(,) k_CtR (x,y ) is Var(k^CtR(,))=12Re[|k^C(,)|2]+⋯ ( k_CtR (x,y ) )= 12Re\! \E [ | k_C (x,y ) |^2 ]+·s ⋯+[(k^C(,))2]−2([k^C(,)])2. ·s+E [ ( k_C (x,y ) )^2 ]-2 (E [ k_C (x,y ) ] )^2 \. (9) For the variance analysis, we need to compute [|k^C(,)|2]E [ | k_C (x,y ) |^2 ] and [(k^C(,))2]E [ ( k_C (x,y ) )^2 ]. By applying Lemma 1 and Lemma 4, we obtain the following bounds: [|k^C(,)|2]≤⟨,⟩2p+2D((⟨,⟩2+ [ | k_C (x,y ) |^2 ]≤ ,y ^2p+ 2D ( ( ,y ^2+ ⋯+∥22∥22−∑i=1dxi2yi2)p−⟨,⟩2p), \ ·s+ \|x \|_2^2 \|y \|_2^2- _i=1^dx_i^2y_i^2 )^p- ,y ^2p ), (10) and [(k^C(,))2]≤⟨,⟩2p+2D((2⟨,⟩2− [ ( k_C (x,y ) )^2 ]≤ ,y ^2p+ 2D ( (2 ,y ^2- ⋯−∑i=1dxi2yi2)p−⟨,⟩2p). ·s- _i=1^dx_i^2y_i^2 )^p- ,y ^2p ). (11) Substituting Equations (10), (11), and (8) into Equation (9) yields Var(k^CtR(,))≤2p+1−2D‖22p‖22p. ( k_CtR (x,y ) )≤ 2^p+1-2D\, \|x \|_2^2p \|y \|_2^2p. (12) Time Complexity. To compute the sketch according to Equation (7), each complex CountSketch can be computed in O(nnz())O(nnz (x )) time. The convolution of the p sketches is implemented using FFT in O(pDlogD)O(pD D) time. Therefore, the overall time complexity of Complex-to-Real (CtR) TensorSketch applied to ⊗px p is O(p(nnz()+DlogD))O\! (p(nnz (x )+D D) ). ∎ Remark 3. Theorem 2, exhibits an improved exponential dependence on p in the variance bound of the estimate, decreasing from the 3p/D3^p/D bound for real TensorSketch [pham2013fast] to 2p/D2^p/D, while retaining the same input-sparsity sketching time. A tabular comparison among the baselines on the variance bound and running time is presented in Table 1. Bounds derived in the Lemma 4 play a key role in the proof of Theorem 2. The sketch described in this lemma is a complex-valued variant of the classical AMS sketch ALON1999137 for polynomial kernel approximation 10.5555/1347082.1347163, braverman2010ams4wiseindependenceproduct. In particular, the variance analysis of CtR TensorSketch reduces to analyzing products of independently randomized linear sketches. Each bucket of CtR TensorSketch behaves like a product of p independent AMS-type sketches of the form Zsj()=∑i=1dxisj(i)Z_s_j(x)= _i=1^dx_is_j(i). When expanding the moments [|Z|2]E[|Z|^2] and [Z2]E[Z^2] in the proof of Theorem 2, the resulting expressions factor across the p independent random functions drawn uniformly from fourth roots of unity. Lemma 4 provides an exact evaluation of moments under fourth roots-of-unity random functions, allowing the full degree-p moment to be written as the p-th power of a single-sketch expression. Lemma 4. Let ,∈ℝdx,y ^d, let p≥1p≥ 1 be an integer, and let s1,…,sp:[d]→1,ω,ω2,ω3s_1,…,s_p:[d]→\1,ω,ω^2,ω^3\ be independent random functions, each taking values uniformly from the four fourth roots of unity (as in Lemma 1). Define Z=∏j=1pZsj()Zsj()¯, Z\;=\; _j=1^pZ_s_j(x)\, Z_s_j(y), where Zsj() Z_s_j(x) =∑i=1dxisj(i), = _i=1^dx_i\,s_j(i), Zsj() Z_s_j(y) =∑i=1dyisj(i). = _i=1^dy_i\,s_j(i). Then, [Z] [Z ] =⟨,⟩p, = ,y ^p, Var(Z) (Z ) ≤2p‖22p‖22p. ≤ 2^p\, \|x \|_2^2p\, \|y \|_2^2p. Proof. Let Z=∏j=1pZsj()Zsj()¯,Zsj()=∑i=1dxisj(i).Z= _j=1^pZ_s_j(x)\, Z_s_j(y), Z_s_j(x)= _i=1^dx_is_j(i). Expectation. For a fixed j, [Zsj()Zsj()¯]=∑i,kxiyk[sj(i)sj(k)¯].E [Z_s_j(x) Z_s_j(y) ]= _i,kx_iy_k\,E [s_j(i) s_j(k) ]. By Lemma 1, [sj(i)sj(k)¯]=0E [s_j(i) s_j(k) ]=0 for i≠ki≠ k and equals 11 for i=ki=k. Hence, [Zsj()Zsj()¯]=⟨,⟩.E [Z_s_j(x) Z_s_j(y) ]= ,y . Independence across j then gives [Z]=∏j=1p⟨,⟩=⟨,⟩p.E [Z ]= _j=1^p ,y = ,y ^p. Second moments. Since the sign functions are independent across j, [|Z|2] [ |Z |^2 ] =∏j=1p[|Zsj()Zsj()¯|2], = _j=1^pE [|Z_s_j(x) Z_s_j(y)|^2 ], [Z2] [Z^2 ] =∏j=1p[(Zsj()Zsj()¯)2]. = _j=1^pE [ (Z_s_j(x) Z_s_j(y) )^2 ]. Expanding one factor and using the fourth–moment structure of fourth roots-of-unity random functions (Lemma 1), only index configurations that form pairs survive. The detailed calculation yields moments as follows [|Zsj()Zsj()¯|2]=⟨,⟩2+‖22‖22−∑i=1dxi2yi2,E [|Z_s_j(x) Z_s_j(y)|^2 ]= ,y ^2+ \|x \|_2^2 \|y \|_2^2- _i=1^dx_i^2y_i^2, and [(Zsj()Zsj()¯)2]=2⟨,⟩2−∑i=1dxi2yi2.E [ (Z_s_j(x) Z_s_j(y) )^2 ]=2 ,y ^2- _i=1^dx_i^2y_i^2. Therefore, taking p-times product of the terms give moment bounds for Z [|Z|2] [|Z|^2] =(⟨,⟩2+‖22‖22−∑i=1dxi2yi2)p, = ( ,y ^2+ \|x \|_2^2 \|y \|_2^2- _i=1^dx_i^2y_i^2 )^p, [Z2] [Z^2 ] =(2⟨,⟩2−∑i=1dxi2yi2)p. = (2 ,y ^2- _i=1^dx_i^2y_i^2 )^p. Variance bound. Using Equation 4, and applying the Cauchy-Schwarz inequality to upper bound ⟨,⟩2≤‖22‖22 ,y ^2≤ \|x \|_2^2 \|y \|_2^2, together with the fact that ∑ixi2yi2≥0 _ix_i^2y_i^2≥ 0, we obtain Var(Z)≤2p‖22p‖22p.Var (Z )≤ 2^p\, \|x \|_2^2p \|y \|_2^2p. ∎ 5 Experiments Experimental setup. Figure 1: KL divergence on the COD–RNA dataset. We report KL divergence between the exact degree-p polynomial kernel and the kernel reconstructed from sketch features for p∈3,5,7p∈\3,5,7\. Methods include Real and CtR (complex-to-real) Gaussian and Rademacher JL sketches, as well as Real and CtR TensorSketch. Results are averaged over 20 independent trials. Sketch dimension is varied as D∈d,3d,5dD∈\d,3d,5d\. Figure 2: Wall-clock sketch construction time on the COD–RNA dataset. We compare Real and CtR (complex-to-real) Gaussian and Rademacher JL sketches, with Real and CtR TensorSketch for p∈3,7,10p∈\3,7,10\. The sketch dimension is varied as D∈d,3d,5dD∈\d,3d,5d\. Each point reports the average sketch construction time over 2020 independent trials, measured on identical normalized input data. Here as well all dense JL-type sketches exhibit nearly overlapping construction times across sketch dimensions D and polynomial degrees p. We compare our proposed CtR TensorSketch against standard TensorSketch [pham2013fast]. Along with this we included JL-type polynomial sketches (Gaussian and Rademacher [kar2012random]) together with their CtR variants [pmlr-v206-wacker23a] for comparison. All methods are implemented by us directly from the original algorithmic descriptions given in the previous sections. All experiments were run on Ubuntu 22.04.4 using an Intel® Core™ i9-14900K processor (24 cores, 32 threads) with 32 GB RAM. Reported runtimes are wall-clock times for sketch construction only. Baselines. The sketching methods compared in this section are listed below, along with brief descriptions. • Real Gaussian/Rademacher Sketch [kar2012random]: Dense JL-type polynomial sketches using i.i.d. Gaussian or Rademacher projection matrices. • Real TensorSketch [pham2013fast, pham2025tensorsketchfastscalable]: Standard CountSketch-based TensorSketch for polynomial kernels with real random sign functions. • CtR Gaussian/Rademacher Sketch [pmlr-v206-wacker23a]: complex-to-real JL-type sketches using complex Gaussian or complex Rademacher projections followed by CtR conversion, as proposed by [pmlr-v206-wacker23a]. • CtR TensorSketch (Our Proposal): Our CtR TensorSketch obtained by replacing real random signs with random function drawn uniformly from the four fourth roots of unity and applying complex-to-real conversion within the corresponding constructions. Datasets. We use both synthetic and real-world datasets. Synthetic data. For the synthetic experiments, we generate n=3000n=3000 random vectors in ℝdR^d with d=2d=2. Each coordinate is drawn independently from a standard Gaussian distribution, and all vectors are ℓ2 _2–normalized. We consider polynomial kernels of degrees p∈10,15,20p∈\10,15,20\ and vary the sketch dimension as D∈100,300,500D∈\100,300,500\. Real–world data. We evaluate all methods on two real-world datasets: MAGIC Gamma Telescope [magic_gamma_telescope_159] (d=10d=10 real-valued features) and COD-RNA [articlecodrna] (d=8d=8 numerical attributes). All inputs are treated as real-valued vectors and ℓ2 _2-normalized prior to sketching. For each dataset, we subsample up to n=3000n=3000 points (or use the full dataset if smaller). Each configuration (dataset, degree p, sketch dimension D) is averaged over 2020 independent random trials. Comparison Metrics. We evaluate sketch quality and efficiency using two metrics. • KL divergence. To assess how well a sketch preserves the structure of the degree-p polynomial kernel matrix, we compute the KL divergence between the exact kernel K and its sketch-based approximation K K. Since polynomial kernels are nonnegative, we normalize both matrices into discrete distributions Pij=Kij∑a,bKab,Qij=K^ij∑a,bK^ab,P_ij= K_ij _a,bK_ab, Q_ij= K_ij _a,b K_ab, and report KL(K∥K^)=∑i,jPijlog(PijQij+ε),KL(K\,\|\, K)= _i,jP_ij \! ( P_ijQ_ij+ ), where ε>0 >0 ensures numerical stability. Lower values indicate better kernel preservation. • Wall–clock time. We measure efficiency by the wall–clock time required to construct sketch features for a given method and sketch dimension D, averaged over multiple independent trials. Results. Figures 1 and 2 compare approximation error (KL divergence) and sketch construction time on the COD–RNA dataset across polynomial degrees and sketch dimensions. We evaluate the proposed CtR TensorSketch against the standard TensorSketch, along with real and CtR Gaussian and Rademacher JL-type sketches. CtR TensorSketch consistently achieves the lowest KL divergence across all degrees and sketch dimensions, indicating more accurate kernel approximation than real TensorSketch and dense JL-type methods. This behavior matches our theoretical guarantee of improved variance established earlier. In terms of runtime, CtR TensorSketch retains the input-sparsity cost O(p(nnz()+DlogD))O(p(nnz (x )+D D)), comparable to standard TensorSketch. In contrast, JL-type methods incur higher dense projection cost O(pDd)O(pDd). Additional experimental results are provided in Section C of the Appendix. 6 Conclusion We introduced a complex-to-real (CtR) variant of TensorSketch that provides a compression algorithm for high-dimensional polynomial-kernel and tensor datasets. To the best of our knowledge, CtR constructions were only known for dense JL-type TensorSketch [pham2013fast] due to pmlr-v206-wacker23a. Our results demonstrate that Complex-to-Real (CtR) variants achieve improved variance bounds compared to their real counterparts [pham2013fast], matching those obtained in pmlr-v206-wacker23a. Moreover, these variants preserve the input-sparsity running time of the original method, making it a preferred choice over pmlr-v206-wacker23a, for high-dimensional sparse data. Two directions remain open for further investigation. First, it is unclear whether employing higher-order roots of unity can yield stronger higher-moment concentration, leading to tighter (ϵ,δ)(ε,δ) -approximation guarantees for the sketch. Second, it remains to be explored whether the advantages of Complex-to-Real (CtR) constructions can be extended to other kernel families and randomized feature maps beyond polynomial kernels. Together, these questions point to a broader potential for complex-valued sketch design. Acknowledgements This research was partially supported by the ANRF under Project No. ANRF/ECRG/2024/001063/ENS. We gratefully acknowledge ANRF for this support. References Improving TensorSketch Using Complex Random Variables (Supplementary Material) Appendix A Proofs of Theorems in Section 4 See 1 Proof. Since s(i)s(i) is uniformly distributed over 1,ω,ω2,ω3=1,i,−1,−i\1,ω,ω^2,ω^3\=\1,i,-1,-i\, we have [s(i)]=1+i−1−i4=0. [s(i) ]= 1+i-1-i4=0. All four values have unit magnitude, and therefore [|s(i)|2]=|1|2+|i|2+|−1|2+|−i|24=1. [ |s(i) |^2]= |1 |^2+ |i |^2+ |-1 |^2+ |-i |^24=1. For the second complex moment, [s(i)2]=12+i2+(−1)2+(−i)24=1−1+1−14=0. [s(i)^2 ]= 1^2+i^2+(-1)^2+(-i)^24= 1-1+1-14=0. Finally, for any i≠ji≠ j, independence of s(i)s(i) and s(j)s(j) implies [s(i)s(j)¯]=[s(i)][s(j)¯]=0. [s(i)\, s(j)]=E[s(i)]\,E [ s(j) ]=0. ∎ A.1 Unbiasedness and Variance Analysis of CtR TensorSketch We establish the theoretical guarantees of CtR TensorSketch. Specifically, we show that the estimator is unbiased for the inner product estimation and subsequently provide a variance analysis. See 2 Proof. According to Definition 5, we have k^CtR(,):=Rek^C(,) k_CtR (x,y ):=Re \ k_C (x,y ) \, and therefore we first analyze k^C(,) k_C (x,y ). Recall that k^C(,):=ΦC(⊗p)ΦC(⊗p)¯⊤ k_C (x,y ):= _C(x p) _C(y p) , where ΦC(⊗p):=⊗p _C(x p):=Cx p and ΦC(⊗p):=⊗p _C(y p):=Cy p according to Definition 8. In particular, ⊗p,⊗p∈ℂD/2Cx p,Cy p ^D/2 denote the complex CountSketch of the vectors ⊗p,⊗p∈ℝdpx p,y p ^d^p, constructed using the derived functions H:[d]p↦[D/2]H:[d]^p [D/2] and S:[d]p→1,ω,ω2,ω3S:[d]^p→\1,ω,ω^2,ω^3\ defined as H(i1,…,ip) H(i_1,…,i_p) =(∑j=1phj(ij))modD/2, = ( _j=1^ph_j(i_j) ) D/2, S(i1,…,ip) S(i_1,…,i_p) =∏j=1psj(ij). = _j=1^ps_j(i_j). In the proof, we denote X:=⊗pX:=x p and Y:=⊗pY:=y p. Let u,v∈[d]pu,v∈[d]^p be the multi-indices corresponding to the entries of X,Y∈ℝdpX,Y ^d^p. For a multi-index u=(j1,…,jp)u=(j_1,…,j_p), the entry of X is defined as Xu=∏k=1pxjkX_u= _k=1^px_j_k, where the associated linearized index is u=1+∑k=1p(jk−1)∏k′=1p−kdu=1+ _k=1^p(j_k-1) _k =1^p-kd. The entries of Y are defined analogously. We begin by analyzing k^C(,) k_C (x,y ) as defined below, k^C(,) k_C (x,y ) :=ΦC(X)ΦC(Y)¯⊤=⟨X,Y¯⟩, := _C(X) _C(Y) = , CY , =∑u,v∈[d]pXuYvS(u)S(v)¯ 1[H(u)=H(v)], = _u,v∈[d]^pX_u\,Y_v\,S(u)\, S(v)\ 1 [H(u)=H(v) ], =⟨X,Y⟩+∑u≠vXuYvS(u)S(v)¯ 1[H(u)=H(v)]. =\ X,Y \!+\! _u≠ vX_u\,Y_v\,S(u)\, S(v)\,1 [H(u)=H(v) ]. Using Lemma 1, we have [S(u)S(v)¯]=0E [S(u)\, S(v) ]=0 for all u≠vu≠ v. Hence, [k^C(,)]=⟨X,Y⟩=⟨⊗p,⊗p⟩=⟨,⟩p. [ k_C (x,y ) ]= X,Y = p,y p = ,y ^p. According to Definition 5, we know k^CtR(,):=Rek^C(,) k_CtR (x,y ):=Re \ k_C (x,y ) \, we have [k^CtR(,)]=[Rek^C(,)]=⟨,⟩p. [ k_CtR (x,y ) ]=E [Re \ k_C (x,y ) \ ]= ,y ^p. (13) Using Equation 5, the variance of k^CtR(,) k_CtR (x,y ) is Var(k^CtR(,))=12Re[|k^C(,)|2]+[(k^C(,))2]−2([k^C(,)])2. ( k_CtR (x,y ) )= 12Re\! \E [ | k_C (x,y ) |^2 ]+E [ ( k_C (x,y ) )^2 ]-2 (E [ k_C (x,y ) ] )^2 \. (14) For the variance analysis, we need to compute [|k^C(,)|2]E [ | k_C (x,y ) |^2 ] and [(k^C(,))2]E [ ( k_C (x,y ) )^2 ]. Therefore by expanding |k^C(,)|2| k_C (x,y )|^2 we get: |k^C(,)|2:=|⟨C⊗p,C⊗p¯⟩|=⟨C⊗p,C⊗p¯⟩⟨C⊗p¯,C⊗p⟩, | k_C (x,y )|^2:=| p, Cy p |= p, Cy p Cx p,Cy p , (15) = = (⟨X,Y⟩+∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])(⟨X,Y⟩+∑u≠vXuYvS(u)¯S(v)[H(u)=H(v)]), ( X,Y + _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] )\!\! ( X,Y + _u≠ vX_uY_v S(u)S(v)1 [H(u)=H(v) ] ), (16) = = ⟨X,Y⟩2+⟨X,Y⟩(∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)]+∑u≠vXuYvS(u)¯S(v)[H(u)=H(v)])⋯ X,Y ^2+ X,Y ( _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ]+ _u≠ vX_uY_v S(u)S(v)1 [H(u)=H(v) ] )·s +|(∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])|2. + | ( _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] ) |^2. (17) Therefore, by applying expectation [|(∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])|2] \! [ | ( _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] ) |^2 ] =[∑u1≠v1u2≠v2Xu1Yv1Xu2Yv2S(u1)S(v1)¯S(u2)¯S(v2)[H(u1)=H(v1)][H(u2)=H(v2)]], =E\! [ _ subarraycu_1≠ v_1\\ u_2≠ v_2 subarrayX_u_1Y_v_1X_u_2Y_v_2S(u_1) S(v_1) S(u_2)S(v_2)1 [H(u_1)=H(v_1) ]1 [H(u_2)=H(v_2) ] ], (18) =∑u1≠v1u2≠v2[Xu1Yv1Xu2Yv2S(u1)S(v1)¯S(u2)¯S(v2)]⋅[[H(u1)=H(v1)][H(u2)=H(v2)]], = _ subarraycu_1≠ v_1\\ u_2≠ v_2 subarrayE\! [X_u_1Y_v_1X_u_2Y_v_2S(u_1) S(v_1) S(u_2)S(v_2) ]·E[1 [H(u_1)=H(v_1) ]1 [H(u_2)=H(v_2) ]], (19) ≤2D∑u1≠v1u2≠v2[Xu1Yv1Xu2Yv2S(u1)S(v1)¯S(u2)¯S(v2)], ≤ 2D _ subarraycu_1≠ v_1\\ u_2≠ v_2 subarrayE\! [X_u_1Y_v_1X_u_2Y_v_2S(u_1) S(v_1) S(u_2)S(v_2) ], (20) ≤2D∑u1≠v1[|Xu1||Yv1||Xu2||Yv2|S(u1)S(v1)¯S(u2)¯S(v2)], ≤ 2D _u_1≠ v_1E\! [|X_u_1|\,|Y_v_1|\,|X_u_2|\,|Y_v_2|\,S(u_1) S(v_1) S(u_2)S(v_2) ], (21) =2D[|(∑u≠v∈[d]p|Xu||Yv|S(u)S(v)¯)|2]. = 2DE\! [ | ( _u≠ v∈[d]^p|X_u|\,|Y_v|\,S(u) S(v) ) |^2 ]. (22) We now bound the above expression using the second-moment bound from Lemma 4, given in Equation (66). For completeness, we first restate the bound, [|(∑u,v∈[d]p|Xu||Yv|S(u)S(v)¯)|2]=(⟨,⟩2+‖22‖22−∑i=1dxi2yi2)p. [ | ( _u,v∈[d]^p|X_u|\,|Y_v|\,S(u) S(v) ) |^2 ]= ( ,y ^2+\|x\|_2^2\,\|y\|_2^2- _i=1^dx_i^2y_i^2 )^p. (23) Now, we expand the term |(∑u,v∈[d]p|Xu||Yv|S(u)S(v)¯)|2 | ( _u,v∈[d]^p|X_u|\,|Y_v|\,S(u) S(v) ) |^2 from the above Equation (23) as follows, |(∑u,v∈[d]p|Xu||Yv|S(u)S(v)¯)|2=|∑u∈[d]p|Xu||Yv|+∑u,v∈[d]pu≠v|Xu||Yv|S(u)S(v)¯|2, | ( _u,v∈[d]^p|X_u|\,|Y_v|\,S(u) S(v) ) |^2= | _u∈[d]^p|X_u|\,|Y_v|+ _ subarraycu,v∈[d]^p\\ u≠ v subarray|X_u|\,|Y_v|\,S(u) S(v) |^2, (24) =(∑u∈[d]p|Xu||Yu|+∑u,v∈[d]pu≠v|Xu||Yv|S(u)S(v)¯)(∑w∈[d]p|Xw||Yw|+∑w,z∈[d]pw≠z|Xw||Yz|S(w)S(z)¯)¯, = ( _u∈[d]^p|X_u|\,|Y_u|\,\;+\; _ subarraycu,v∈[d]^p\\ u≠ v subarray|X_u|\,|Y_v|\,S(u) S(v) ) ( _w∈[d]^p|X_w|\,|Y_w|\,\;+\; _ subarraycw,z∈[d]^p\\ w≠ z subarray|X_w|\,|Y_z|\,S(w) S(z) ), (25) =(∑u∈[d]p|Xu||Yu|+∑u,v∈[d]pu≠v|Xu||Yv|S(u)S(v)¯)(∑w∈[d]p|Xw||Yw|+∑w,z∈[d]pw≠z|Xw||Yz|S(w)¯S(z)), = ( _u∈[d]^p|X_u|\,|Y_u|\,\;+\; _ subarraycu,v∈[d]^p\\ u≠ v subarray|X_u|\,|Y_v|\,S(u) S(v) ) ( _w∈[d]^p|X_w|\,|Y_w|\,\;+\; _ subarraycw,z∈[d]^p\\ w≠ z subarray|X_w|\,|Y_z|\, S(w)S(z) ), (26) =∑u,w∈[d]p|Xu||Yu||Xw||Yw|+∑u∈[d]p∑w,z∈[d]pw≠z|Xu||Yu||Xw||Yz|S(w)¯S(z)+⋯ = _u,w∈[d]^p|X_u|\,|Y_u|\,|X_w|\,|Y_w|\,+ _u∈[d]^p _ subarraycw,z∈[d]^p\\ w≠ z subarray|X_u|\,|Y_u|\,|X_w|\,|Y_z|\, S(w)S(z)+·s ⋯+∑u,v∈[d]pu≠v∑w∈[d]p|Xu||Yv||Xw||Yw|S(u)S(v)¯+∑u,v∈[d]pu≠v∑w,z∈[d]pw≠z|Xu||Yv||Xw||Yz|S(u)S(v)¯S(w)¯S(z). ·s+ _ subarraycu,v∈[d]^p\\ u≠ v subarray _w∈[d]^p|X_u|\,|Y_v|\,|X_w|\,|Y_w|\,S(u) S(v)\,+ _ subarraycu,v∈[d]^p\\ u≠ v subarray _ subarraycw,z∈[d]^p\\ w≠ z subarray|X_u|\,|Y_v|\,|X_w|\,|Y_z|\,S(u) S(v)\, S(w)S(z). (27) Using Lemma 1, we have [S(u)S(v)¯]=0E [S(u)\, S(v) ]=0 for all u≠vu≠ v. Hence, ∑u,v∈[d]pu≠v∑w∈[d]p|Xu||Yv||Xw||Yw|[S(u)S(v)¯]=0, _ subarraycu,v∈[d]^p\\ u≠ v subarray _w∈[d]^p|X_u|\,|Y_v|\,|X_w|\,|Y_w|\,E\! [S(u) S(v) ]=0, (28) Similarly, for all w≠zw≠ z, ∑u∈[d]p∑w,z∈[d]pw≠z|Xu||Yu||Xw||Yz|[S(w)¯S(z)]=0, _u∈[d]^p _ subarraycw,z∈[d]^p\\ w≠ z subarray|X_u|\,|Y_u|\,|X_w|\,|Y_z|\,E\! [ S(w)S(z) ]=0, (29) Substitute this in Equation (23), we get [|∑u,v∈[d]p|Xu||Yv|S(u)S(v)¯|2] [ | _u,v∈[d]^p|X_u|\,|Y_v|\,S(u) S(v) |^2 ] =∑u,w∈[d]p|Xu||Yu||Xw||Yw|+⋯ = _u,w∈[d]^p|X_u|\,|Y_u|\,|X_w|\,|Y_w|\ +·s ⋯+[∑u,v∈[d]pu≠v∑w,z∈[d]pw≠z|Xu||Yv||Xw||Yz|S(u)S(v)¯S(w)¯S(z)], ·s+E [ _ subarraycu,v∈[d]^p\\ u≠ v subarray _ subarraycw,z∈[d]^p\\ w≠ z subarray|X_u|\,|Y_v|\,|X_w|\,|Y_z|\,S(u) S(v)\, S(w)S(z) ], (30) =⟨X,Y⟩2+|(∑u≠v|Xu||Yv|S(u)S(v)¯)|2. = X,Y ^2+E | ( _u≠ v|X_u|\,|Y_v|\,S(u) S(v) ) |^2. (31) We conclude that, [|∑u≠v|Xu||Yv|S(u)S(v)¯|2] [ | _u≠ v|X_u|\,|Y_v|\,S(u) S(v) |^2 ] =[|∑u,v∈[d]p|Xu||Yv|S(u)S(v)¯|2]−⟨,⟩2p. =E [ | _u,v∈[d]^p|X_u|\,|Y_v|\,S(u) S(v) |^2 ]- ,y ^2p. (32) Putting value of [|(∑u,v∈[d]p|Xu||Yv|S(u)S(v)¯)|2]E [ | ( _u,v∈[d]^p|X_u|\,|Y_v|\,S(u) S(v) ) |^2 ], we get [|∑u≠v|Xu||Yv|S(u)S(v)¯|2] [ | _u≠ v|X_u|\,|Y_v|\,S(u) S(v) |^2 ] =(⟨,⟩2+‖22‖22−∑i=1dxi2yi2)p−⟨,⟩2p. = ( ,y ^2+\|x\|_2^2\,\|y\|_2^2- _i=1^dx_i^2y_i^2 )^p- ,y ^2p. (33) Put this value in Equation 22, [|∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)]|2]≤2D((⟨,⟩2+‖22‖22−∑i=1dxi2yi2)p−⟨,⟩2p). \! [ | _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] |^2 ]≤ 2D ( ( ,y ^2+\|x\|_2^2\,\|y\|_2^2- _i=1^dx_i^2y_i^2 )^p- ,y ^2p ). (34) Now we have the second moment as follows, [|k^C(,)|2] [| k_C (x,y )|^2 ] =⟨X,Y⟩2+[|(∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])|2], = X,Y ^2+E\! [ | ( _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] ) |^2 ], (35) ≤⟨,⟩2p+2D((⟨,⟩2+‖22‖22−∑i=1dxi2yi2)p−⟨,⟩2p). ≤ ,y ^2p+ 2D ( ( ,y ^2+\|x\|_2^2\,\|y\|_2^2- _i=1^dx_i^2y_i^2 )^p- ,y ^2p ). (36) Next, we analyze [(k^C(,))2]E\! [( k_C (x,y ))^2 ]. Expanding (k^C(,))2( k_C (x,y ))^2 gives: (k^C(,))2:=(⟨C⊗p,C⊗p¯⟩)2=⟨C⊗p,C⊗p¯⟩⟨C⊗p,C⊗p¯⟩, ( k_C (x,y ))^2:= ( p, Cy p )^2= p, Cy p p, Cy p , (37) =(⟨X,Y⟩+∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])(⟨X,Y⟩+∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)]), = ( X,Y + _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] ) ( X,Y + _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] ), (38) =⟨X,Y⟩2+2⟨X,Y⟩(∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])+(∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])2, = X,Y ^2+2 X,Y ( _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] )+ ( _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] )^2, (39) Therefore, by applying expectation [(∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])2] \! [ ( _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] )^2 ] (40) =[∑u1≠v1u2≠v2Xu1Yv1Xu2Yv2S(u1)S(v1)¯S(u2)S(v2)¯[H(u1)=H(v1)][H(u2)=H(v2)]], =E\! [ _ subarraycu_1≠ v_1\\ u_2≠ v_2 subarrayX_u_1Y_v_1X_u_2Y_v_2S(u_1) S(v_1)S(u_2) S(v_2)1 [H(u_1)=H(v_1) ]1 [H(u_2)=H(v_2) ] ], (41) =∑u1≠v1u2≠v2[Xu1Yv1Xu2Yv2S(u1)S(v1)¯S(u2)S(v2)¯]⋅[[H(u1)=H(v1)][H(u2)=H(v2)]], = _ subarraycu_1≠ v_1\\ u_2≠ v_2 subarrayE\! [X_u_1Y_v_1X_u_2Y_v_2S(u_1) S(v_1)S(u_2) S(v_2) ]·E[1 [H(u_1)=H(v_1) ]1 [H(u_2)=H(v_2) ]], (42) ≤2D[(∑u≠v∈[d]pXuYvS(u)S(v)¯)2]. ≤ 2DE\! [ ( _u≠ v∈[d]^pX_u\,Y_v\,S(u) S(v) )^2 ]. (43) By similar analysis and using second-moment bound of Lemma 4 given in Equation (75), we get [(∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])2]≤2D((2⟨,⟩2−∑i=1dxi2yi2)p−⟨,⟩2p). \! [ ( _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] )^2 ]≤ 2D ( (2 ,y ^2- _i=1^dx_i^2y_i^2 )^p- ,y ^2p ). (44) Therefore, we have, [(k^C(,))2] [ ( k_C (x,y ) )^2 ] =⟨X,Y⟩2+[(∑u≠vXuYvS(u)S(v)¯[H(u)=H(v)])2], = X,Y ^2+E\! [ ( _u≠ vX_uY_vS(u) S(v)1 [H(u)=H(v) ] )^2 ], (45) ≤⟨,⟩2p+2D((2⟨,⟩2−∑i=1dxi2yi2)p−⟨,⟩2p). ≤ ,y ^2p+ 2D ( (2 ,y ^2- _i=1^dx_i^2y_i^2 )^p- ,y ^2p ). (46) Therefore for final variance we put Equation (36) and Equation (46) in Equation (14),we get Var(k^CtR(,)) ( k_CtR (x,y ) ) ≤12[⟨,⟩2p+2D((⟨,⟩2+∥22∥22−∑i=1dxi2yi2)p−⟨,⟩2p) ≤ 12 [ ,y ^2p+ 2D ( ( ,y ^2+\|x\|_2^2\,\|y\|_2^2- _i=1^dx_i^2y_i^2 )^p- ,y ^2p ) +⟨,⟩2p+2D((2⟨,⟩2−∑i=1dxi2yi2)p−⟨,⟩2p)−2⟨,⟩2p], + ,y ^2p+ 2D ( (2 ,y ^2- _i=1^dx_i^2y_i^2 )^p- ,y ^2p )-2 ,y ^2p ], (47) we can upper bound the above equation by using inequality ⟨,⟩2≤‖22‖22 ,y ^2≤\|x\|^2_2\|y\|^2_2, then ≤12[2D((2‖22‖22)p−‖22p‖22p)+2D((2‖22‖22)p−‖22p‖22p)], ≤ 12 [ 2D ( (2\|x\|_2^2\,\|y\|_2^2 )^p-\|x\|_2^2p\,\|y\|_2^2p )+ 2D ( (2\|x\|_2^2\,\|y\|_2^2 )^p-\|x\|_2^2p\,\|y\|_2^2p ) ], (48) ≤12[2p+1−2D‖22p‖22p+2p+1−2D‖22p‖22p], ≤ 12 [ 2^p+1-2D\|x\|_2^2p\,\|y\|_2^2p+ 2^p+1-2D\|x\|_2^2p\,\|y\|_2^2p ], (49) ≤2p+1−2D‖22p‖22p. ≤ 2^p+1-2D\|x\|_2^2p\,\|y\|_2^2p. (50) ∎ We begin by stating a lemma that serves as a key step to bound the second moments of the estimator in the proof of Theorem 2 . See 4 Proof. First, we consider the expectation. For each j, we note that [Zsj()Zsj()¯] \! [Z_s_j(x)\, Z_s_j(y) ] =[(∑i=1dxisj(i))(∑k=1dyksj(k)¯)], =E\! [ ( _i=1^dx_i\,s_j(i) ) ( _k=1^dy_k\, s_j(k) ) ], (51) =∑i=1d∑k=1dxiyk[sj(i)sj(k)¯], = _i=1^d _k=1^dx_iy_k\,E[s_j(i) s_j(k)], (52) =∑i=1dxiyi[|sj(i)|2]+∑i≠kxiyk[sj(i)sj(k)¯], = _i=1^dx_iy_i\,E[|s_j(i)|^2]+ _i≠ kx_iy_k\,E[s_j(i) s_j(k)], (53) =⟨,⟩, = ,y , (54) where, [sj(i)sj¯(k)]=0,∀i≠kE[s_j(i) s_j(k)]=0,∀ i≠ k and [|sj(i)|2]=1,∀i∈[d]E[|s_j(i)|^2]=1,∀\ i∈[d] as given in Lemma 1. Since the functions sjs_j are independent across different j, we have [Z]=∏j=1p[Zsj()Zsj()¯]=⟨,⟩p. [Z]= _j=1^pE[Z_s_j(x) Z_s_j(y)]= ,y ^p. (55) Next, to bound the variance, Var(Z)=12Re[|Z|2]+[(Z)2]−2|[Z]|2. (Z )= 12Re\! \E[|Z|^2]+E[(Z)^2]-2|E[Z]|^2 \. (56) Because the random hash functions sjs_j are independent across different j, we may write [|Z|2]=∏j=1p[|(Zsj()Zsj()¯)|2] and [(Z)2]=∏j=1p[(Zsj()Zsj()¯)2] [|Z|^2]= _j=1^pE\! [| (Z_s_j(x)\, Z_s_j(y) )|^2 ] and E[(Z)^2]= _j=1^pE\! [ (Z_s_j(x)\, Z_s_j(y) )^2 ] (57) For each j, expanding the square gives [|(Zsj()Zsj()¯)|2] \! [| (Z_s_j(x)\, Z_s_j(y) )|^2 ] =[(∑i=1dxisj(i))(∑k=1dyksj(k)¯)(∑i=1dxisj(i)¯)(∑k=1dyksj(k))], =E\! [ ( _i=1^dx_i\,s_j(i) ) ( _k=1^dy_k\, s_j(k) ) ( _i=1^dx_i\, s_j(i) ) ( _k=1^dy_k\,s_j(k) ) ], (58) =∑i=1d∑i′=1d∑k=1d∑k′=1dxixi′ykyk′[sj(i)sj(k)¯sj(i′)¯sj(k′)]. = _i=1^d _i =1^d _k=1^d _k =1^dx_i\,x_i \,y_k\,y_k \;E\! [s_j(i) s_j(k) s_j(i )s_j(k ) ]. (59) Observing that [sj(i)sj(k)¯sj(i′)¯sj(k′)]E[s_j(i) s_j(k) s_j(i )s_j(k )] is nonzero only when the indices form pairs (including the possibility that all four are identical), we have [sj(i)sj(k)¯sj(i′)¯sj(k′)]=1,if i=k=i′=k′,1,if i=k≠i′=k′,1,if i=i′≠k=k′,0,otherwise. [s_j(i) s_j(k) s_j(i )s_j(k )]= cases1,&if i=k=i =k ,\\[4.0pt] 1,&if i=k≠ i =k ,\\[4.0pt] 1,&if i=i ≠ k=k ,\\[4.0pt] 0,&otherwise. cases (60) The contribution from terms with i=k=i′=k′i=k=i =k is ∑i=1dxi2yi2. _i=1^dx_i^2y_i^2. (61) Terms with i=k≠i′=k′i=k≠ i =k contribute ∑i≠i′xiyixi′yi′=(∑i=1dxiyi)2−∑i=1dxi2yi2=⟨,⟩2−∑i=1dxi2yi2. _i≠ i x_iy_i\,x_i y_i = ( _i=1^dx_iy_i )^2- _i=1^dx_i^2y_i^2= ,y ^2- _i=1^dx_i^2y_i^2. (62) Finally, for i=i′≠k=k′i=i ≠ k=k we obtain ∑i≠kxi2yk2=‖22‖22−∑i=1dxi2yi2. _i≠ kx_i^2y_k^2=\|x\|_2^2\,\|y\|_2^2- _i=1^dx_i^2y_i^2. (63) Thus, summing these contributions, we have [|(Zsj()Zsj())|2] \! [| (Z_s_j(x)\,Z_s_j(y) )|^2 ] =∑i=1dxi2yi2+(⟨,⟩2+‖22‖22−2∑i=1dxi2yi2), = _i=1^dx_i^2y_i^2+ ( ,y ^2+\|x\|_2^2\,\|y\|_2^2-2 _i=1^dx_i^2y_i^2 ), (64) =⟨,⟩2+‖22‖22−∑i=1dxi2yi2. = ,y ^2+\|x\|_2^2\,\|y\|_2^2- _i=1^dx_i^2y_i^2. (65) Substituting this bound into Equation (57) yields [|Z|2]=(⟨,⟩2+‖22‖22−∑i=1dxi2yi2)p, [|Z|^2]= ( ,y ^2+\|x\|_2^2\,\|y\|_2^2- _i=1^dx_i^2y_i^2 )^p, (66) For each j, expanding the square gives [(Zsj()Zsj()¯)2] \! [ (Z_s_j(x)\, Z_s_j(y) )^2 ] =[(∑i=1dxisj(i))(∑k=1dyksj(k)¯)(∑i=1dxisj(i))(∑k=1dyksj(k)¯)], =E\! [ ( _i=1^dx_i\,s_j(i) ) ( _k=1^dy_k\, s_j(k) ) ( _i=1^dx_i\,s_j(i) ) ( _k=1^dy_k\, s_j(k) ) ], (67) =∑i=1d∑i′=1d∑k=1d∑k′=1dxixi′ykyk′[sj(i)sj(k)¯sj(i′)sj(k′)¯]. = _i=1^d _i =1^d _k=1^d _k =1^dx_i\,x_i \,y_k\,y_k \;E\! [s_j(i) s_j(k)s_j(i ) s_j(k ) ]. (68) Observing that [sj(i)sj(k)¯sj(i′)sj(k′)¯]E[s_j(i) s_j(k)s_j(i ) s_j(k )] is nonzero only when the indices form pairs (including the possibility that all four are identical), we have [sj(i)sj(k)¯sj(i′)sj(k′)¯]=1,if i=k=i′=k′,1,if i=k≠i′=k′,1,if i=k′≠i′=k,0,otherwise. [s_j(i) s_j(k)s_j(i ) s_j(k )]= cases1,&if i=k=i =k ,\\[4.0pt] 1,&if i=k≠ i =k ,\\[4.0pt] 1,&if i=k ≠ i =k,\\[4.0pt] 0,&otherwise. cases (69) The contribution from terms with i=k=i′=k′i=k=i =k is ∑i=1dxi2yi2. _i=1^dx_i^2y_i^2. (70) Terms with i=k≠i′=k′i=k≠ i =k contribute ∑i≠i′xiyixi′yi′=(∑i=1dxiyi)2−∑i=1dxi2yi2=⟨,⟩2−∑i=1dxi2yi2. _i≠ i x_iy_i\,x_i y_i = ( _i=1^dx_iy_i )^2- _i=1^dx_i^2y_i^2= ,y ^2- _i=1^dx_i^2y_i^2. (71) Finally, for i=k′≠i′=ki=k ≠ i =k we obtain ∑i≠i′xiyixi′yi′=(∑i=1dxiyi)2−∑i=1dxi2yi2=⟨,⟩2−∑i=1dxi2yi2. _i≠ i x_iy_i\,x_i y_i = ( _i=1^dx_iy_i )^2- _i=1^dx_i^2y_i^2= ,y ^2- _i=1^dx_i^2y_i^2. (72) Thus, summing these contributions, we have [(Zsj()Zsj())2] \! [ (Z_s_j(x)\,Z_s_j(y) )^2 ] =∑i=1dxi2yi2+(2⟨,⟩2−2∑i=1dxi2yi2), = _i=1^dx_i^2y_i^2+ (2 ,y ^2-2 _i=1^dx_i^2y_i^2 ), (73) =2⟨,⟩2−∑i=1dxi2yi2. =2 ,y ^2- _i=1^dx_i^2y_i^2. (74) Substituting this bound into Equation (57) yields [(Z)2]=(2⟨,⟩2−∑i=1dxi2yi2)p, [(Z)^2]= (2 ,y ^2- _i=1^dx_i^2y_i^2 )^p, (75) which completes the proof since Var(Z) (Z ) =12Re[|Z|2]+[(Z)2]−2|[Z]|2, = 12Re\! \E[|Z|^2]+E[(Z)^2]-2|E[Z]|^2 \, (76) =12(⟨,⟩2+‖22‖22−∑i=1dxi2yi2)p+(2⟨,⟩2−∑i=1dxi2yi2)p−2⟨,⟩2p. = 12 \ ( ,y ^2+\|x\|_2^2\,\|y\|_2^2- _i=1^dx_i^2y_i^2 )^p+ (2 ,y ^2- _i=1^dx_i^2y_i^2 )^p-2 ,y ^2p \. (77) Using the Cauchy–Schwarz inequality, ⟨,⟩2≤‖22‖22 ,y ^2≤\|x\|_2^2\,\|y\|_2^2, and noting that ∑i=1dxi2yi2≥0 _i=1^dx_i^2y_i^2≥ 0, it follows that Var(Z)≤2p‖22p‖22p. (Z )≤ 2^p\|x\|_2^2p\|y\|_2^2p. (78) ∎ Appendix B proof of CountSketch Estimator with Fourth Roots of Unity Theorem 5. Let Φ denote the CountSketch with the hash function taking values in the fourth roots of unity. For any vectors ,∈ℝdx,y ^d, the corresponding sketch vector is defined by Φ()j=∑i=1dδjiσ(i)xi,Φ()j=∑i=1dδjiσ(i)yi,∀j∈[D/2]. (x)_j= _i=1^d _jiσ(i)x_i\,, (y)_j= _i=1^d _jiσ(i)y_i, ∀ j∈[D/2]. (79) where, • h:[d]→[D/2]h:[d]→[D/2] assigning each coordinate independently and uniformly to one of D/2D/2 buckets, and • σ:[d]→1,ω,ω2,ω3σ:[d]→\1,ω,ω^2,ω^3\ are drawn independently and uniformly at random, • δji=1,if h(i)=j,0,otherwise. _ji= cases1,&if h(i)=j,\\ 0,&otherwise. cases then the estimator k^(,):=ReΦ()⊤Φ()¯, k(x,y):=Re \ (x) (y) \, satisfies [k^(,)] \! [ k(x,y) ] =⟨,⟩, = ,y , (80) Var(k^(,)) \! ( k(x,y) ) =1D(‖22‖22+⟨,⟩2−2∑i=1dxi2yi2). = 1D\! (\|x\|_2^2\|y\|_2^2+ ,y ^2-2 _i=1^dx_i^2y_i^2 ). (81) Proof. Since σ(i)σ(i) is uniformly distributed over 1,ω,ω2,ω3=1,i,−1,−i\1,ω,ω^2,ω^3\=\1,i,-1,-i\ and satisfies [σ(i)]=0,[|σ(i)|2]=1,[(σ(i))2]=0E[σ(i)]=0,\,E\! [\,|σ(i)|^2\, ]=1,\,E\! [(σ(i))^2 ]=0. Consider now the complex inner product of the sketch: k^C(,) k_C(x,y) :=Φ()⊤Φ()¯=∑j=1D/2(∑i:h(i)=jσ(i)xi)(∑m:h(m)=jσ(m)¯ym) := (x) (y)= _j=1^D/2 ( _i:\,h(i)=jσ(i)x_i ) ( _m:\,h(m)=j σ(m)\,y_m ) (82) =∑i=1d|σ(i)|2xiyi+∑i≠mh(i)=h(m)σ(i)σ(m)¯xiym. = _i=1^d|σ(i)|^2x_iy_i+ _i≠ m1\h(i)=h(m)\\,σ(i) σ(m)\,x_iy_m. (83) Taking expectation with respect to h and σ, we get [k^C(,)]=∑i=1d[|σ(i)|2]xiyi=∑i=1dxiyi=⟨,⟩. \! [ k_C(x,y) ]= _i=1^dE\! [|σ(i)|^2 ]x_iy_i= _i=1^dx_iy_i= ,y . (84) Since the right-hand side of (84) is real, the same unbiasedness holds for the estimator, which takes the real part of the complex sketch: [k^(,)] \! [ k(x,y) ] =[Rek^C(,)]=Re[k^C(,)]=⟨,⟩. =E [Re \ k_C(x,y) \ ]=Re \E [ k_C(x,y) ] \= ,y . (85) To find the exact variance of the estimator, we leverage the Complex-to-Real (CtR) framework stated in our paper. Now, we work out with [|k^C(,)|2]E\! [| k_C(x,y)|^2 ]. [|k^C(,)|2] \! [| k_C(x,y)|^2 ] =[|Φ()⊤Φ()¯|2]=[|∑j=1D/2(∑i=1dδjiσ(i)xi)(∑m=1dδjmσ(m)¯ym)|2], =E\! [| (x) (y)|^2 ]=E\! [ | _j=1^D/2 ( _i=1^d _jiσ(i)x_i ) ( _m=1^d _jm σ(m)\,y_m ) |^2 ], (86) =∑j,j′D/2∑i,m,i′,m′d[δjiδjmδj′iδj′m′][σ(i)σ(m)¯σ(i′)¯σ(m′)]xiymxi′ym′. = _j,j ^D/2 _i,m,i ,m ^dE\! [ _ji _jm _j i _j m ]\;E\! [σ(i) σ(m)\, σ(i )σ(m ) ]\;x_iy_mx_i y_m . (87) Hence, [σ(i)σ(m)¯σ(i′)¯σ(m′)]≠0⇔i,m,i′,m′appear in pairs,E\! [σ(i) σ(m)\, σ(i )σ(m ) ]≠ 0 \i,m,i ,m \\ appear in pairs, 1. With surviving index configurations, when =′j=j : • ==′=′:[δji][|σ(i)|4]xi2yi2i=m=i =m :E[ _ji]E[|σ(i)|^4]x_i^2y_i^2, • =≠′=′:[δjiδji′][|σ(i)|2|σ(i′)|2]xiyixi′yi′i=m≠ i =m :E[ _ji _ji ]E[|σ(i)|^2|σ(i )|^2]x_iy_ix_i y_i , • =′≠=′:[δjiδjm][|σ(i)|2|σ(m)|2]xi2ym2i=i ≠ m=m :E[ _ji _jm]E[|σ(i)|^2|σ(m)|^2]x_i^2y_m^2, 2. With surviving index configurations, when ≠′j≠ j : • =≠′=′:[δjiδj′i′][|σ(i)|2|σ(i′)|2]xiyixi′yi′i=m≠ i =m :E[ _ji _j i ]E[|σ(i)|^2|σ(i )|^2]x_iy_ix_i y_i , Thus, [|k^C(,)|2] \! [| k_C(x,y)|^2 ] =∑i=1dxi2yi2+2D∑i≠i′(xiyixi′yi′+xi2yi′2)+D−2D∑i≠i′xiyixi′yi′, = _i=1^dx_i^2y_i^2+ 2D _i≠ i (x_iy_ix_i y_i +x_i^2y_i ^2)+ D-2D _i≠ i x_iy_ix_i y_i , (88) =⟨,⟩2+2D∑i≠i′xi2yi′2. = ,y ^2+ 2D _i≠ i x_i^2y_i ^2. (89) Now, similarly we get [(k^C(,))2] \! [( k_C(x,y))^2 ] =∑i=1dxi2yi2+4D∑i≠i′(xiyixi′yi′)+D−2D∑i≠i′xiyixi′yi′, = _i=1^dx_i^2y_i^2+ 4D _i≠ i (x_iy_ix_i y_i )+ D-2D _i≠ i x_iy_ix_i y_i , (90) =⟨,⟩2+2D∑i≠i′(xiyixi′yi′). = ,y ^2+ 2D _i≠ i (x_iy_ix_i y_i ). (91) Now, substitute values of [|k^C(,)|2]E\! [| k_C(x,y)|^2 ] and [(k^C(,))2]E\! [( k_C(x,y))^2 ] and simplifying, we get Var(k^(,)) \! ( k(x,y) ) =1D(‖22‖22+⟨,⟩2−2∑i=1dxi2yi2). = 1D (\|x\|_2^2\|y\|_2^2+ ,y ^2-2 _i=1^dx_i^2y_i^2 ). (92) Hence, from the above variance bound, it follows that CountSketch with complex random variables does not provide any variance reduction over its real-valued counterpart. ∎ Appendix C Further Experiments C.1 Extended Evaluation on Variance and Time in different setups In this section, we present additional experimental results that complement the main empirical evaluation reported in Section 5. These experiments extend the comparison between real and complex-to-real (CtR) sketching constructions across additional datasets, higher polynomial degrees, and varied embedding dimension. We report both approximation quality (via KL divergence between exact and sketched kernels) and wall-clock sketch construction time. Together, these results provide a more detailed view of the variance and construction time trade-offs of CtR TensorSketch relative to its real-valued counterparts and JL-type baselines. Figure 3: KL divergence on the MAGIC Gamma Telescope dataset. We compare Real and complex-to-real (CtR) Gaussian and Rademacher JL sketches, together with Real and CtR TensorSketch. Results are shown for polynomial degrees p∈3,5,7p∈\3,5,7\ and sketch dimensions D∈d,3d,5dD∈\d,3d,5d\ with n=3000n=3000 standardized and ℓ2 _2-normalized samples. Bars report the mean KL divergence over 20 independent trials. Figure 4: Wall-clock sketch construction time on the MAGIC dataset. Methods compared include Real and complex-to-real (CtR) Gaussian and Rademacher JL sketches, with Real and CtR TensorSketch. Results are shown for polynomial degrees p∈3,5,7p∈\3,5,7\ and sketch dimensions D∈d,3d,5dD∈\d,3d,5d\ with n=3000n=3000 standardized and ℓ2 _2-normalized samples. Each point reports the mean runtime over 20 independent trials for feature sketch construction only. Figure 5: KL divergence on the synthetic dataset (d=2d=2). We compare Real and complex-to-real (CtR) Gaussian and Rademacher JL sketches, together with Real and CtR TensorSketch, for approximating degree-p polynomial kernels on synthetic Gaussian data (n=3000n=3000, dimension d=2d=2, standardized and ℓ2 _2-normalized). Results are shown for polynomial degrees p∈10,15,20p∈\10,15,20\ and sketch dimensions D∈100,300,500D∈\100,300,500\. Bars report the mean KL divergence between the exact kernel and the sketch-based approximation over 20 independent trials. Figure 6: Wall-clock sketch construction time on a synthetic dataset. Methods compared include Real and complex-to-real (CtR) Gaussian and Rademacher JL sketches, together with Real and CtR TensorSketch. Results are shown for polynomial degrees p∈10,15,20p∈\10,15,20\ and sketch dimensions D∈100,300,500D∈\100,300,500\ with n=3000n=3000 standardized and ℓ2 _2-normalized samples in dimension d=2d=2. Each point reports the mean runtime over 20 independent trials for feature sketch construction only (log-scale on the y-axis). C.2 Evaluation on Frobenius Normalized Relative Error We also test our proposed method along with other baselines using the Frobenius normalized relative error, which is defined as ‖K−K^‖F/‖K‖F\|K- K\|_F/\|K\|_F. In this formula, K represents the exact kernel matrix, while K K represents the approximated kernel matrix produced by the sketching method. Frobenius error metric is a standard benchmark in kernel approximation and is widely used to measure sketching accuracy wacker2022improved, pmlr-v206-wacker23a. The datasets and experimental setups used here are identical to our previous experiments that evaluated accuracy using the KL divergence metric (Section 5 and Appendix C.1). Testing with this alternative metric ensures that our method performs well irrespective of the choice of error metric. Degree (p) Method Variance 15 Real Gaussian 4.41×1044.41× 10^4 Real Rademacher 3.20×1013.20× 10^1 CtR Gaussian 6.35×1016.35× 10^1 CtR Rademacher 4.16×1034.16× 10^3 Real TensorSketch 2.84×10−52.84× 10^-5 CtR TensorSketch (Ours) 1.57 × 10-5 20 Real Gaussian 3.80×1013.80× 10^1 Real Rademacher 6.78×1036.78× 10^3 CtR Gaussian 2.76×1012.76× 10^1 CtR Rademacher 1.58×1011.58× 10^1 Real TensorSketch 1.28×10−21.28× 10^-2 CtR TensorSketch (Ours) 5.05 × 10-5 25 Real Gaussian 2.65×1002.65× 10^0 Real Rademacher 9.06×1039.06× 10^3 CtR Gaussian 1.82×1011.82× 10^1 CtR Rademacher 2.58×1012.58× 10^1 Real TensorSketch 3.68×10−13.68× 10^-1 CtR TensorSketch (Ours) 3.58 × 10-4 30 Real Gaussian 7.24×1007.24× 10^0 Real Rademacher 1.00×1001.00× 10^0 CtR Gaussian 5.44×1025.44× 10^2 CtR Rademacher 3.21×1033.21× 10^3 Real TensorSketch 1.27×10−31.27× 10^-3 CtR TensorSketch (Ours) 4.59 × 10-6 Table 2: Variance of the Frobenius normalized relative error (‖K−K^‖F/‖K‖F\|K- K\|_F/\|K\|_F) evaluated on the real-world MAGIC dataset (n=1000,d=10n=1000,d=10) with a compressed sketch dimension of D=128D=128. Results are aggregated over 20 independent trials across high-degree polynomial kernels (p∈15,20,25,30p∈\15,20,25,30\). CtR TensorSketch strictly outperforms Real TensorSketch in estimation stability, confirming the variance reduction. C.3 Evaluation on Downstream Tasks To demonstrate practical application beyond kernel matrix approximation, we evaluate our method on downstream binary classification tasks using a linear Support Vector Machine (SVM). The datasets used in these experiments are generated using standard ML libraries (such as scikit-learn) to create two highly overlapping classes with non-linear decision boundaries. We test two separate configurations under extreme compression to a sketch dimension of D=64D=64. We compare the sketching methods against the Exact Polynomial Kernel, which computes the full kernel matrix using the mathematical formula K(x,y)=(x⊤y)pK(x,y)=(x y)^p without any compression or approximation. We evaluate the techniques using two metrics: test classification accuracy and total execution time. The time metric tracks the combined end-to-end wall-clock time required for both data compression and the subsequent classifier training. We report this total time to fully reflect the complete workload required to obtain the final classification model. As shown in Table 3 and Table 4, all approximation methods experience a drop in accuracy compared to the exact kernel due to the tight bottleneck of the compressed dimension. However, our proposed CtR TensorSketch consistently achieves the highest accuracy among all baselines while requiring the shortest total execution time. Method Accuracy Time (s) Exact Poly Kernel 0.8756 16.0724 TensorSketch (Real) 0.4911 0.2390 TensorSketch (CtR) [Ours] 0.5289 0.1479 JL (CtR Rademacher) 0.5011 3.8279 JL (CtR Gaussian) 0.4978 3.7708 Table 3: Downstream Classification Performance (Linear SVM) for Setup 1. Evaluated on a synthetic dataset (n=3000n=3000) with dimension d=20d=20. We approximate a polynomial kernel of degree p=15p=15 using a sketch dimension of D=64D=64. Our proposed method achieves the best sketching accuracy in the shortest time. Method Accuracy Time (s) Exact Poly Kernel 0.8567 8.5344 TensorSketch (Real) 0.4978 0.2763 TensorSketch (CtR) [Ours] 0.5044 0.1733 JL (CtR Rademacher) 0.4567 3.3204 JL (CtR Gaussian) 0.4889 3.2954 Table 4: Downstream Classification Performance (Linear SVM) for Setup 2. Evaluated on a synthetic dataset (n=3000n=3000) with dimension d=10d=10. We approximate a polynomial kernel of higher degree p=25p=25 using a sketch dimension of D=64D=64. C.4 Variance Analysis through Numerical Table The primary purpose of this section is to provide a clear and direct validation of our experimental findings from Section 5. In our main evaluation, we rely on visual plots to demonstrate the accuracy and computational time of our proposed CtR TensorSketch. Graphs, especially those that use a logarithmic scale, can visually compress the performance gap between methods. Looking at the exact numbers, we can clearly confirm our previous experiments and highlight the advantage of our method. Specifically, these numbers reveal how our approach successfully reduces the variance from 3p/D3^p/D in traditional methods to 2p/D2^p/D in our complex-to-real design. As shown in Table 5, our proposed CtR TensorSketch consistently achieves lower variance compared to all the baselines across all degrees. Degree (p) Method Variance 15 Real Gaussian 0.2634 Real Rademacher 0.3839 CtR Gaussian 0.2192 CtR Rademacher 0.2029 Real TensorSketch 0.5496 CtR TensorSketch (Ours) 0.1739 20 Real Gaussian 0.3545 Real Rademacher 0.4673 CtR Gaussian 0.6996 CtR Rademacher 0.5246 Real TensorSketch 0.3932 CtR TensorSketch (Ours) 0.3888 25 Real Gaussian 0.3158 Real Rademacher 0.2560 CtR Gaussian 0.4395 CtR Rademacher 0.2294 Real TensorSketch 0.3209 CtR TensorSketch (Ours) 0.2367 30 Real Gaussian 0.4785 Real Rademacher 0.2932 CtR Gaussian 0.3266 CtR Rademacher 0.2883 Real TensorSketch 0.2443 CtR TensorSketch (Ours) 0.1939 Table 5: Variance of KL Divergence for high-degree polynomial kernels. Evaluated on a positive orthant synthetic dataset (n=1000,d=10n=1000,d=10) with extreme compression (D=64D=64), aggregated over 20 independent trials. CtR TensorSketch strictly dominates Real TensorSketch in estimation variance across all polynomial degrees, empirically validating the theoretical 2p/D2^p/D scaling advantage.