Locality Sensitive Hashing in Hyperbolic Space

Locality Sensitive Hashing (LSH) is a hashing scheme in which similar items are more likely to hash to the same bucket. Indyk and Motwani, in their seminal work [27, 25], developed LSH for general metric spaces and provided hash function constructions for the Hamming distance. The LSH framework provides an efficient data structure for the Approximate Nearest Neighbor Search (ANNS) problem, which is particularly challenging in the high-dimensional regimes. Since then, LSH has been developed for various similarity measures such as ${\mathrm{\ell }}_1$ [1, 35], ${\mathrm{\ell }}_2$ [18, 2], ${\mathrm{\ell }}_p$ [18], cosine similarity [15], Jaccard similarity [11, 10], angular (spherical) distances [47], etc. With support for similarity search, data independence, mild dependence on data dimension, and simplicity in implementation, LSH serves as an algorithmic tool for classical problems such as clustering, near-duplicate detection, and efficient information retrieval for very large datasets. LSH has also been used widely in database and machine learning pipelines in industry [22, 12, 16, 44, 42, 17, 32, 4, 7]. For a nice survey, please refer to [2].

Hyperbolic geometry is a non-Euclidean geometry. It drops the parallel postulate of Euclidean geometry. In hyperbolic geometry, for any line and a point not on that line, there are at least two distinct lines through the point that do not intersect the original line. Hyperbolic geometry has received increasing attention in modern data analysis because it naturally models hierarchical, tree-like, and exponentially branching structures that arise throughout science and technology. Unlike Euclidean spaces, where volume grows polynomially, hyperbolic spaces exhibit exponential volume growth, matching the combinatorial growth of hierarchies, taxonomies, knowledge graphs, and high-dimensional networks with latent tree-like structure. This geometric compatibility yields embeddings with dramatically reduced distortion: objects that are far apart in a hierarchy can be faithfully separated, and geodesics efficiently capture ancestry, flow, and community structure. As a result, hyperbolic embeddings now underpin state-of-the-art methods in network representation learning, hierarchical clustering, and discrete geometric data analysis [38, 24, 21, 48, 19, 13, 14, 49, 30, 46]. It is also the most important space representing the negative curvature phenomena arising from complex geometric situations [34, 45]. Namely, according to Gromov’s hyperbolization program [23], spaces with complicated topology tend to resemble spaces with negative curvature. Hyperbolic geometry provides the fundamental model for negatively curved spaces, offering the canonical framework in which their curvature and geometric behavior are understood.

Moreover, the hyperbolic plane offers an unusually tractable analytic structure: its constant negative curvature provides explicit formulas for distance, geodesics, convexity, and isometries, enabling efficient optimization for machine-learning tasks. Its rich isometry group and uniformization properties are directly linked to complex analysis and conformal representations, providing tools that scale from theoretical guarantees to practical algorithms for dimensionality reduction, visualization, and manifold learning.

These applications provide strong motivation for developing efficient algorithms in hyperbolic space, and in this paper, we study locality-sensitive hashing in hyperbolic space, a topic that, to our best knowledge, has not been studied before.

In this section, we show the LSH construction in ${ℍ}^2$ that yields Theorem 4. We adopt the Poincaré disk model for points in ${ℍ}^2$, defined on the unit disk . The space is equipped with the metric $ds=\frac{2|dz|}{1-{|z|}^2}$, which induces the following distance function:

Geometrically, geodesics in $𝔻$ correspond to circular arcs orthogonal to the unit circle boundary $C_{\mathrm{\infty }}$, with diameters passing through the origin as special cases.

The locality sensitive hash functions we use are random geodesics in the Poincaré disk model. For this purpose we require a well-defined probability measure on the set of all geodesics to ensure that random sampling of geodesics is meaningful and invariant under hyperbolic isometries. Luckily, the classical integral geometry [43] provides the needed theory. Integral geometry studies the measures on the spaces of all lines, or totally geodesic planes in Euclidean, hyperbolic, or spherical spaces, such that the measures are invariant under isometries. The fundamental theorems of Poincaré, Cauchy and Crofton show that one can recover basic geometric quantities such as lengths, areas, or volumes of an object by averaging over the set of geodesics intersecting the object.

In the hyperbolic plane, the isometrically invariant probability measure on the space of all geodesics is given by the kinematic measure [43]. Take a point $q$ and a non-zero tangent vector $v$ at $q$ in the hyperbolic plane. The polar coordinate of a geodesic $\gamma$ with respect to $(q,v)$ is a pair $(t,\theta )$ where $t$ is the hyperbolic distance from $\gamma$ to $q$ and $\theta$ is the angle of the perpendicular from $q$ to $\gamma$ with the tangent vector $v$. The kinematic measure in geodesic polar coordinates $(t,\theta )$ is

The remarkable feature of the above formula is that it is independent of the choices of $q$ and $v$. This measure reflects the fact that the volume element in hyperbolic geometry grows exponentially with $t$.

Since the space of all geodesics in ${ℍ}^2$ has infinite total measure, we limit our choices of geodesics in the family of all geodesics intersecting a hyperbolic disk $B(0,R)$ of radius $R$ centered at the origin. By Crofton’s formula, the measure of the set of geodesics that intersect a convex body is proportional to the perimeter of that body. Specifically, the corresponding measure is obtained by integrating the invariant element $d\mu$ over the disk:

Now, we choose $R$ such that the disk $B(0,R)$ centered at the origin with radius $R$ contains all points of $P$. Our final bound on the performance parameter $\rho$ does not depend on $R$, and in practice, one can choose $R$ to be a sufficiently large value. Note that in practice, hyperbolic embedding typically avoids using points that are far from the origin due to practical concerns on resolution and numerical instability [36].

A geodesic $\gamma$ is chosen with probability density function $\frac{\cosh (t)}{2\pi \sinh (R)}$ to be at polar coordinate $(t,\theta )$ where $t$ is taken from $[0,R]$ and $\theta$ is taken from $[0,2\pi )$. In the Poincaré disk model, geodesics appear as Euclidean circular arcs, so $\gamma$ is a circular arc with Euclidean center $c\in {ℝ}^2$ and radius $b>0$. Specifically, define $c=\coth (t)e^{i\theta }$ (here convert $e^{i\theta }$ from the complex plane to the Euclidean plane) and $b=1/\sinh t$. Define a hash function ${ℍ}^2\to \{-1,1\}$ by

where $x$ is the Euclidean coordinate of a point in ${ℍ}^2$. We can also rewrite the geodesic $\gamma$ using the Minkowski model with inner product $⟨z,z⟩=-z_1^2+z_2^2+z_3^2$. Let the geodesic $\gamma$ be given by $\{x\in {ℍ}^2|⟨x,u⟩=0\}$ for a point $u$ in the de Sitter space $\{y\in {ℝ}^{2,1}|⟨y,y⟩=1\}$:

Note that $(t,\theta )$ is polar coordinate with respect to the point $q=(1,0,0)$ and $v=(0,1,0)$. Then the formula for $h_{\gamma }(x)$ is $\text{sgn}⟨x,u⟩.$ This is exactly the same formula as in the Euclidean space. The advantage of this formula over the above one is that it is intrinsic to hyperbolic geometry, i.e., invariant under hyperbolic isometry.

With this setup, we present the hash function in ${ℍ}^2$ below.

For brevity, we write $h(x)$ for $h_{\gamma }(x)$ when it is clear from context. Figure 1 demonstrates how the hash function works as a separator.

The key to proving Theorem 4 is a characterization of the collision probability under our hashing scheme, given in the following lemma.

In this section, we prove the dimension reduction result (Theorem 12). We use the Poincaré half-space model of the hyperbolic space ${ℍ}^d$. Each point in ${ℍ}^d$ is represented by a pair $(z,x)$ with $z\in {ℝ}^{+}$ and $x\in {ℝ}^{d-1}$. The hyperbolic distance between two points $p_1=(z_1,x_1)$ and $p_2=(z_2,x_2)$ is defined by the following function

where $\Vert x_1-x_2\Vert$ is the Euclidean distance between $x_1,x_2$. For a set of $n$ points $P=\{p_i=(z_i,x_i)\}\subseteq {ℍ}^d$, we would like to perform dimension reduction to generate $P^{\prime }=\{p_i^{\prime }\}$ in ${ℍ}^k$ and keep the pairwise distances similar. The method by Benjamini and Makarychev [9] is to perform dimension reduction for the Euclidean point set $X=\{x_i\}\subseteq {ℝ}^{d-1}$ by a function $f$ to $X^{\prime }=\{x_i^{\prime }\}\subseteq {ℝ}^{k-1}$ and define $P^{\prime }=\{p_i^{\prime }\}$ with $p_i^{\prime }=(z_i,x_i^{\prime })$ in ${ℍ}^k$. Below, we present a refined analysis that yields slightly tighter bounds on the distortion. Essentially, we show that the bound on the distortion of the Euclidean distances of $X^{\prime }$ (compared to $X$) carries over to bound the distortion of the hyperbolic distance of $P^{\prime }$ compared to $P$.

The proof follows the main idea in [9], but we have slightly improved analysis and bounds. Specifically, the main idea is to analyze the behavior of the function

which is an increasing function of $r$, when $z_1,z_2>0$. For completeness, we include the full proof in the full version [20].

A distribution $D$ is called $p$-stable [51] with $p\ge 0$ if for any $n$ real numbers $v_1,v_2,\mathrm{\cdots },v_n$ and random independent variables $X_1,X_2,\mathrm{\cdots }{X}_n$ with distribution $D$, ${\sum }_i{v}_i{X}_i$ has the same distribution as the variable ${({\sum }_i{|v_i|}^p)}^{1/p}X$, where $X$ is a random variable with distribution $D$. A Gaussian distribution $𝒩(0,1)$ is $2$-stable. This property has been used to build Euclidean LSH [18], to sketch high-dimensional vectors [26], and for many other applications.

Consider $n$ points $P=\{p_i\}$ in ${ℍ}^d$. We use random projection to project these points to $P^{\prime }=\{p_i^{\prime }\}$ in ${ℍ}^2$ using the method in Section 3. In particular, we take a vector $\overrightarrow{a}$ with dimension $d-1$ where each entry is independently taken from a Gaussian distribution $𝒩(0,1)$. For $p_i=(z_i,x_i)$ in the upper half space model, we have $p_i=(z_i,x_i^{\prime })$ in ${ℍ}^2$, where $x_i^{\prime }=\overrightarrow{a}\cdot x_i$. We then use the LSH mechanism in Section 2 to map these points into buckets. In this section, we analyze the performance of this mechanism.

In our case, the value $\overrightarrow{a}\cdot (p-q)$ for any two vectors $p,q\in {ℝ}^{d-1}$ is also a Gaussian distribution with zero mean and variance ${\Vert p-q\Vert }^2$, or, the distribution of $\Vert p-q\Vert Z$ with $Z\sim 𝒩(0,1)$. Take $f(z)$ to be the probability density function of the absolute value of $𝒩(0,1)$:

To calculate the probability that $p_1,p_2$ are mapped to the same bucket, we take $s=\Vert x_1-x_2\Vert$ and $r=d_{ℍ}(p_1,p_2)=F_{z_1,z_2}(s)$, where

is an increasing function of $s$. Suppose $p_1,p_2$ are mapped to points $p_1^{\prime },p_2^{\prime }$ with $s^{\prime }=\Vert x_1^{\prime }-x_2^{\prime }\Vert$ and $r^{\prime }=d_{ℍ}(p_1^{\prime },p_2^{\prime })=F_{z_1,z_2}(s^{\prime })$. $s^{\prime }$ follows the distribution of $sZ$ with $Z\sim 𝒩(0,1)$. Recall that using the LSH for ${ℍ}^2$ the probability that $p_1^{\prime },p_2^{\prime }$ map to the same bucket is $1-\frac{r^{\prime }}{w}$ where $w=2\sinh R$. Now we have,

We can get an upper bound of $p(r)$.

The lower bound on the performance parameter $\rho$ for LSH of ${\mathrm{\ell }}_p$ metrics has been studied in [37, 40]. Theorem 5.2 in [40] stated that for points of the Hamming cube ${\{0,1\}}^d$ and some , an $(r,R,p_1,p_2)$-sensitive LSH of Hamming distance, with $r=\delta d/c$ and $R=\delta d$, must have

This lower bound can be used to generate a lower bound of $1/c^p$ for ${\mathrm{\ell }}_p$ distance – by using ${\mathrm{\ell }}_p$ distance on the same point set. In particular, it shows that LSH for Euclidean distances has $\rho \ge 1/c^2$. We will use the same construction to show a lower bound on $\rho$ for hyperbolic distance.

We implement our Locality Sensitive Hashing algorithms in both ${ℍ}^2$ and ${ℍ}^d$. We aim to evaluate the LSH performance in terms of $\rho$ in practice.

We first explain the synthetic data we use. For each $d\in \{2,10,100,1000\}$, we generate a dataset of 1000 points by uniformly random sampling in a hyperbolic sphere of the same radius, $R=\log (199)$. This radius was chosen so that after the points are projected to ${ℍ}^2$ and mapped to the Poincare disk model, the points lie inside a Euclidean circle of radius 0.99. Further, we set $r=0.2$ with varying $c=1.5+k$ for $k=0,1\mathrm{\dots }17$ to compute the value of $p_1$ and $p_2$. We illustrate results on $p_1,p_2$ in Figure 3 and $\rho$ in Figure 4 for different choices of $c$ obtained by averaging over 1000 repetitions.

We see that, as expected, $p_1$ and $p_2$ are close to 1 in all cases. Due to the randomness of the data, the probabilities vary slightly in percentage. Since $\rho$ is related to the complements of these probabilities, varying slightly in proportion when so close to 1 causes the significant differences in $\rho$ values that we are observing. However, in all cases, the $\rho$ graphs are below $1/c$. We note that they are not exactly $1/c$ because we considered all pairs of points with distance below $r$ or above $cr$, not exactly equal to $r$ or $cr$, so $1/c$ should be an upper bound of what we obtain.

Next, we use the same data for $d=2$, but vary $c=1+.1k$ for $k=1,\mathrm{\dots }10$ to see if the LSH still works for small values of $c$ even though our proof only has guarantees for $c\ge 1.59$. From Figure 5 we see fairly similar results, just $p_2$ varying less, but LSH still achieves $\rho$ lower than $\frac{1}{c}$.

We also generated some data for $d=1000$, varying $R$ so that the points lie in a circle of Euclidean radii $.9$ and $.999$. We see that $p_1,p_2$ and $\rho$ change a little, but not significantly and still fall well within the theoretical bounds.

Lastly, we generate similar data for $d=2$, but containing 2000 points. Then, we analyze $p_1$ and $p_2$ for pairs of points that are at the boundary of acceptable distances. Precisely, we choose pairs of points that have distance within $[.9r,r]$ or $[cr,1.1cr]$. This gives us an approximately linear graph for $p_2$ and $\rho$ curve that approximates $\frac{1}{c}$.

The experiment setup with a fixed $r$ and a reasonably large constant $c$ (rather than $1+\epsilon$ for a very small $\epsilon$) is motivated by real world applications such as similarity learning which focuses on differentiating data elements that are very similar or sufficiently far away. Many data modalities and data tasks fall into this category. For example, in image recognition the goal is to recognize two pictures with similar scenes or semantics. But for two pictures on completely different subjects, it is sufficient to know that they are far away and the exact values of dissimilarities are less important. One important take-away message from our experiments is that the empirical value of $\rho$ with the synthetic data set is very small with a mild dependency on dimension and decreasing value of $c$. It is much smaller than our theoretical upper bound – partly because we consider all pairs whose distance are often much smaller than $r$ or much greater than $cr$. Since the storage requirement for approximate nearest neighbor data structure built with LSH grows in terms of $O(n^{1+\rho })$ and the query cost grows in terms of $O(n^{\rho })$. The empirical small value of $\rho$ suggests great potential for practical use.