---
# System prepended metadata

title: Implementing singular value decomposition in zero-knowledge proofs
tags: [Templates, Documentation]

---

---
title: 'Halo2-svd: Singular value decomposition in zero-knowledge proofs'
disqus: hackmd

---
Halo2-svd: Singular value decomposition in zero-knowledge proofs
===

## Table of Contents

[TOC]

## Introduction

We present [Halo2-svd](https://github.com/goforashutosh/halo2-svd), an open source library which provides functionality to verify the singular value decomposition (SVD) of a matrix and to perform basic linear algebraic operations in Halo2. Singular value decomposition is one of the most important tools available in linear algebra to analyze a matrix. Our library brings this powerful tool to zero-knowledge proofs, allowing developers an easy and efficient method for performing the linear algebra required for applications ranging from machine learning to decentralized reputation management.

Using singular value decomposition, every matrix $M$ can be written as the product
$$M = U \Sigma V,$$
where $U$ and $V$ are orthogonal matrices and $\Sigma$ is a diagonal matrix with positive entries. Essentially, this decomposition says that a matrix can be viewed as a rotation of coordinates and a rescaling. The orthogonal matrices determine the rotation and the diagonal matrix $\Sigma$ determines the rescaling. With this interpretation, it can be seen why SVD is an important tool in machine learning- by eliminating the small entries in $\Sigma$, one can *compress* the matrix $M$, while retaining its most essential features. This is essentially the idea behind [Principal Component Analysis](https://en.wikipedia.org/wiki/Principal_component_analysis) (PCA), which is technique for dimensionality reduction in machine learning. The singular vectors (these are the columns of $U$ and the rows of $V$) also satisfy special properties, which have been used in algorithms, such as [Eigentrust](https://nlp.stanford.edu/pubs/eigentrust.pdf).

Zero-knowledge (ZK) proofs [4] are a framework that allow a prover to prove that a statement is true, without revealing any additional information beyond the validity of the statement. One of the most popular implementation of a ZK proof system is [Halo2](https://github.com/zcash/halo2), which provides basic primitives for writing ZK circuits. Recently this has been used by [Axiom](https://www.axiom.xyz) to trustlessly access historic on-chain data of Ethereum. The availability of this data naturally prompts interest in the application of statistical and machine learning algorithms. To facilitate this, an efficient implementation of fundamental linear algebraic operations such as matrix multiplication, eigenvalue calculations, and singular value decompositions within a ZK circuit becomes imperative.

There are two major challenges to overcome in order to implement these operations in ZK circuits:
1. Performing matrix multiplication is a very costly operation: for size $n$ square matrices, it naively requires $O(n^3)$ operations and with Strassen's algorithm $O(n^{2.8})$ [3].
2. One needs to convert floating point numbers to the prime field $\mathbb{F}_p$, since all numbers in a ZK circuit are represented in such a field. This requires implementing some sort of quantization scheme and a corresponding fixed point arithmetic.

For the first point above, the key observation is that one simply needs to *verify* the result of a matrix multiplication inside the ZK circuit, not compute it. This can be done with only $O(n^2)$ operations with very high probability, see **Lemma 1** below. The second obstacle is overcome by using the ZKFixedPointChip [1], which is a library that implements fixed point arithmetic for Halo2. (See also the blog post [1] for an explanation of the Fixed Point Chip and the quantization method that it uses to convert floating point numbers to elements of the prime field.)

The prime field $\mathbb{F}_p$ we will work with is constructed from a $254$-bits prime described in [5]. In the following we will use $M_{n \times m}(\mathbb{R})$ and  $M_{n \times m}(\mathbb{F}_p)$ to denote $n\times m$ matrices valued in $\mathbb{R}$ and $\mathbb{F_p}$ respectively.

## Efficient algorithms for matrix multiplication and SVD

We start by describing the algorithm for the verification of matrix multiplication. Consider three matrices $A \in M_{n\times k}(\mathbb{F}_p)$, $B \in M_{k \times m}(\mathbb{F}_p)$, and $C \in M_{n\times m}(\mathbb{F}_p)$. Suppose that we wish to verify that $C= AB$ inside the ZK-circuit. Instead of naively multiplying the matrices $A$ and $B$ inside the circuit and checking the result against $C$, we use the following variant of Freivalds' randomized algorithm [8] to verify matrix multiplication (interestingly, this variant appears in Sec 2.2 of [9]):
1. Sample a random field element $r \in \mathbb{F}_p$. 
2. Construct the vector $v = (1, r, \dots, r^{m-1})$.
3. Compute the quantities $A(Bv)$ and $Cv$. 
4. If $ABv = Cv$, then the test passes and we have verified that $AB=C$; otherwise, it fails.

The element $r$ is randomly sampled by the Prover using the Fiat-Shamir transform [6]. Note that it is important that the Prover commits to matrices $A$, $B$ and $C$ before the Fiat-Shamir transform is used to sample $r$. The computation of the multiplication result $C$ is performed by the Prover outside the circuit and he can use any algorithm for this purpose. In our code, we simply use the trivial multiplication algorithm.

The most important aspect of the above algorithm is that it only requires $nm + mk + kn$ multiplications (i.e., one multiplication per element of the matrices), compared to $O(nmk)$ multiplications for the trivial matrix multiplication algorithm.

The following lemma implies that if the claimed multiplication result is incorrect, then the test fails with high probability.

**Lemma 1.** Let $A \in M_{n\times k}(\mathbb{F}_p)$, $B \in M_{k \times m}(\mathbb{F}_p)$ and $C \in M_{n\times m}(\mathbb{F}_p)$. Let $r \in \mathbb{F}_p$ be chosen uniformly at random over $\mathbb{F}_p$. Define $v = (1, r, \dots, r^{m-1})$. If $AB \neq C$ (as field matrices), then one has:
$$ \mathbb{P}(Cv = ABv) \leq  \frac{m-1}{p}.$$

*Proof.* This result is a simple consequence of the Schwartz-Zippel lemma [2]. If $AB \neq C$, then there exists $i,j$ such that $c_{i,j} \neq \sum_{l=1}^k a_{i,l}b_{l,j}$. By a direct computation the coefficient $i$ of the vector $Cv - ABv$ is given by:
\begin{align*}
(Cv - ABv)_i = \sum_{s=1}^m r^{s-1} \left( c_{i,s} - \sum_{l=1}^k a_{i, l} b_{l, s} \right). 
\end{align*}
By our assumption $(Cv - ABv)_i$ is a non-zero polynomial in $r$ of degree at most $m-1$. 
Therefore $(Cv - ABv)_i = 0$ is equivalent to a non-zero polynomial being zero at the random field element $r$. The Schwartz-Zippel lemma implies this can happen with probability at most $(m-1)/ p$. Q.E.D.

Note that we have so far only discussed matrix multiplication for matrices over the prime field. What we are really interested in are matrices over real numbers. We will discuss the matrix multiplication for such matrices in the next section. 

Once we have an efficient algorithm for verifying matrix multiplication, we can use it to verify the SVD $(U, \Sigma, V)$ of matrix $M$ by simply checking that $M \approx U\Sigma V$ and that $U^tU \approx Id$ and $V^t V \approx Id$. 

## Design considerations
### Matrix multiplication

Under fixed point precision, there are multiple methods of implementing Freivalds' algorithm [8] for real matrices-- the most straightforward being that we trivially encode all operations using the fixed point chip. For example to check if $C= AB$ (say all matrices with dimension $n \times n$), analogous to the field matrix algorithm described above, we can check that their quantizations satisfy $C_q \approx A_q B_q$ (see next section for a precise definition of the quantization scheme used by the fixed point chip) by iterating the following simple test $k$ times:

1. Sample a random vector $v \in \{-1,1\}^n$ 
2. Compute $u := (C_q - A_q B_q)v$
3. Check if $\Vert u \Vert < \delta$ for some tolerance $\delta$ 

While such an implementation appears to be correct, one needs to be very careful while implementing this using the fixed point chip, especially for the check in Step 3, because the fixed point chip does not check for overflows and a dishonest prover could exploit this fact to cheat. Secondly, it turns out that for a reasonably small probability of error for this algorithm (say $10^{-9}$), the number of iterations required is $k \geq 30$, which is quite significant and leads to a large overhead (between 5-10 times compared to our implementation) in circuit size. Lastly, each iteration would also require a hash to be calculated in circuit, which too is quite expensive. 

Instead, in our library we implement Freivalds' algorithm at a *lower* level- at the level of the field elements. In the fixed point arithmetic of [1], a $P$ precision element $x \in \mathbb{R}$ is represented as the field element $x \cdot 2^P$  in the field $\mathbb{F}_p$. Real number multiplication under this encoding is simply field multiplication followed by a rescaling of the output by $2^P$. Matrix multiplication of real matrices under this encoding can then be viewed as field matrix multiplication followed by a rescaling of each of the elements of the output. Freivalds' algorithm for field matrices can be implemented as described in the previous section. The implementation and analysis of this algorithm is much more efficient and simple than the one for the real numbers, particularly because for field matrices, we do not need to check if the difference vector ($u$ above in Step 2) is *close* to zero, we can simply check that it is zero or not. Further, it is sufficient to use a single random vector for the test in the algorithm, which eliminates the overhead due to multiple iterations. Finally, it turns out that this random vector can also be chosen to be of the form $(1, r, r^2, \cdots, r^{n-1})$ for a random field element $r$ as is shown above in Lemma 1. As a result, we do not require any hashes to be computed in circuit apart from the commitment for the first phase of the ZK-circuit, which generates the element $r$. 

Overall field matrix multiplication can be checked using about ~$9n^2$ cells. As we mentioned above, to obtain the real matrix output from the output of the field matrix multiplication, we need to rescale the matrix. This costs between $60n^2$ to $100n^2$ (depending on the lookup table size and the precision) and turns out to be the much more expensive part of the multiplication procedure. 

### SVD verification

We design the SVD verification procedure such that it does not require rescaling any matrix product, since, as mentioned above, this is a very expensive operation. The bounds, which we choose to check for the verification of the SVD can all be viewed as $\Vert M_1 M_2 - M_3 M_4\Vert_m < \delta$ (see the next section and also [7]). For such bounds, it is sufficient to instead create the unscaled field products of these matrices and simply verify that the difference of their elements is smaller than an appropriately scaled version of $\delta$. In addition to these bounds, we also verify that the elements of $U_q$ and $V_q$ are smaller than $1$ and that the singular values are smaller than the maximum value that can be represented using the fixed point chip. These bounds ensure that the matrix multiplications used do not overflow. 



## The singular value decomposition algorithm

We now describe the full algorithm to check the singular value decomposition of a matrix $M \in M_{n \times m}(\mathbb{R})$ inside the ZK-circuit. The singular value decomposition (SVD) of a matrix $M$ can be written as
$$ M = U \Sigma V,$$
where $U \in O_n(\mathbb{R})$, $V\in O_m(\mathbb{R})$ and $\Sigma \in M_{n\times m}(\mathbb{R})$ is a diagonal matrix with non-negative diagonal entries. Here $O_n(\mathbb{R})$ denotes the set of real orthogonal matrices of size $n \times n$.

There are two sources of numerical error that will arise as one attempts to check such a decomposition within a ZK-circuit. First when one numerically computes the SVD of $M$ using for instance numpy, one will get matrices $\hat U, \hat \Sigma, \hat V$ which are only approximately equal to their exact counterparts $U,\Sigma, V$. Second, one needs to apply a quantization scheme to these matrices in order to be able to represent then in $\mathbb{F}_p$. The quantization of a matrix $X$ is given by:
$$ (X_q)_{i,j} = 2^{-P} \lfloor 2^P X_{i,j} + 1/2\rfloor.$$
In other words, $X_q$ is the matrix with the elements of $X$, up to a precision of $P$ bits. $2^P$ is referred to as the quantization scale. Therefore for any real matrix $X \in M_{n \times m}(\mathbb{R})$, the matrix $2^P X_q$ can be represented as a field matrix in $M_{n \times m}(\mathbb{F}_p)$.


We now describe our complete algorithm to verify the SVD of $M$ in Halo2. We use the observation that given $VV^T = Id$, we have:
$$ M = U \Sigma V \: \:  \Leftrightarrow \: \: M V^{T} = U \Sigma.  $$
Rewriting the SVD condition this way will help us avoid rescaling the matrices while checking the SVD. 

SVD computation and verification in our library follows the following steps. Suppose the prover wishes to compute the SVD of a fixed point quantized matrix $M$ (i.e., $M_q = M$), then the prover:

1. Computes the numerical approximations $\hat U, \hat \Sigma, \hat V$ to the matrices in the singular value decomposition $M = U \Sigma V$ using any library for linear algebra (numpy in the example in our library).
2. Converts the matrices $\hat U, \hat \Sigma, \hat V$ to prime field valued matrices using the quantization procedure described above, namely the matrices:
$$U_{P} := 2^P \hat U_q, \quad  \Sigma_{P} := 2^P \hat \Sigma_q, \quad V_{P} := 2^P \hat V_q. $$
Similarly define $M_{P} := 2^P M_q = 2^P M$ since $M_q= M$.
3. Computes the matrix product $P_1:= M_{P}V^{T}_{P}$ outside the ZK-circuit.
4. Computes the product $P_2:= U_{P} \Sigma_{P}$ inside the ZK-circuit. Note this only requires $n^2$ multiplications since $\Sigma_{P}$ is a diagonal matrix. 
5. Proves the following claims in the ZK-circuit: 
* The matrix multiplication $P_1:= M_{P}V^{T}_{P}$ is correct using the randomized algorithm of the previous section.
* That $||P_1 - P_2||_{m} <2^{2P} \delta_1$. The definition of the norm $||\cdot||_{m}$ and the way $\delta_1$ should be chosen will be explained in the next section.
6. Computes the two products $P_U := U_{P} U_{P}^T$ and $P_V := V_{P} V_{P}^T$ outside the ZK-circuit, then proves inside the ZK-circuit:
* The matrix multiplications $P_U = U_{P} U_{P}^T$ and $P_V = V_{P} V_{P}^T$ are correct.
* That $||2^{2P} \mathrm{Id} - P_U||_{m} < 2^{2P} \delta_2$ and $||2^{2P} \mathrm{Id} - P_V||_{m} < 2^{2P} \delta_2$, where $\delta_2$ is a small error also defined in the next section.

Note that the above algorithm can easily be adapted to check the eigenvalue decomposition of a real symmetric matrix $M$.

## Error analysis of the SVD algorithm.

In this section, we present precise error bounds used in our SVD verification algorithm. For an in-depth error analysis, we refer the reader to [7], where these bounds are proven.

Let us first explain how we define the norm $|| \cdot ||_{m}$ defined over field valued matrices. Given an integer $x$ satisfying $0 \leq x < p$, where $p$ is the prime of our prime field $\mathbb{F}_p$, define the quantity:
$$|x|_{\mathbb{F}_p} := \begin{cases}x & 0 < x \leq \lfloor p/2 \rfloor, \\
p - x & \lfloor p/2 \rfloor \leq x < p. \end{cases}$$

The idea behind this definition is that we want to think of prime field elements past the midpoint of the field as representing negative numbers. For a matrix $X \in M_{m \times n}(\mathbb{F}_p)$, we can now define:
$$ ||X||_{m} := \max_{i,j} |X_{i,j}|_{\mathbb{F}_p}.$$

Now we can state the error analysis results. Let $\epsilon = 2^{-53}$. The results of [7] guarantee that:
* $||P_1 - P_2||_{m} < 2^{P} \max(n,m) (1+ ||M||_2) + 3 \epsilon \cdot 2^{2P} \max(n,m) ||M||_2,$
* $||2^{2P} \mathrm{Id} - P_U||_{m} < 2^{2P}\epsilon \max(n,m) + 3 \cdot 2^{P-1} n.$

So we want to choose $\delta_1 = 2^{-P} \max(n,m) (1+ ||M||_2) + 3   \epsilon \max(n,m) ||M||_2$ and $\delta_2 = \epsilon \max(n,m) + 3 \cdot 2^{-P-1} n$.

To further get a sense of the values of $\delta_1, \delta_2$, let us assume for instance that $M$ is a square matrix of size $n=m=1000$. We choose $P=32$ for the quantization scale and further assume that $||M||_2 < 10^2$. Then one gets roughly that $\delta_1 \approx 2^{-32} \cdot 10^5 \approx 2.3 \cdot 10^{-5}$ and that $\delta_2 \approx 3 \cdot 2^{-33} \cdot 1000 \approx 3.5 \cdot 10^{-7}$. The first error is larger since it is also growing with the matrix norm. 


## Performance of the algorithm

The total number of advice cells and lookup cells used will depend on the parameter $P$ of quantization scale which is denoted by PRECISION_BITS in the code, as well as on the LOOKUP_BITS parameter.

For LOOKUP_BITS =19 and PRECISION_BITS= 32, the number of advice cells grow as $162N^2 = 135N^2 + 27N^2$ (first phase + second phase) with $26N^2$ lookup cells, and for LOOKUP_BITS =19 and PRECISION_BITS= 63 the number of advice cells grows as $228N^2 = 201N^2 + 27N^2$ (first phase + second phase) with $48N^2$ lookup cells. Note here that we use the [Halo2 Challenge API](https://hackmd.io/@axiom/SJw3p-qX3) to implement the Fiat-Shamir transform [6] and sample the randomness for matrix multiplication verification. This API divides the program into two phases. 


## Acknowledgments

This project was conducted as part of the [Axiom Open Source Program](https://www.axiom.xyz/open-source-v2). We are grateful to N. Feng, Y. Sun, and J. Wang for organizing it. It's build on top of ZK circuit primitives provided by [Halo2](https://github.com/privacy-scaling-explorations/halo2.git) and Axiom's [Halo2-base](https://github.com/axiom-crypto/halo2-lib). We also used the fixed-point chip implementation of [1].

## Reference

- [1] ZK Fixed Point Arithmetic with its Application in Machine Learning, https://hackmd.io/JaoTzYgqSY6lXDN4C5alXQ?view
- [2] Schwartz–Zippel Lemma https://en.wikipedia.org/wiki/Schwartz–Zippel_lemma
- [3] Strassen algorithm https://en.wikipedia.org/wiki/Strassen_algorithm
- [4] Zero-knowledge proof https://en.wikipedia.org/wiki/Zero-knowledge_proof
- [5] BN254 For The Rest Of Us https://hackmd.io/@jpw/bn254
- [6] Fiat–Shamir heuristic https://en.wikipedia.org/wiki/Fiat–Shamir_heuristic
- [7] Error bounds for SVD verification in Halo2 https://github.com/goforashutosh/halo2-svd/blob/svd/Error%20Analysis%20for%20SVD.pdf
- [8] Freivalds' algorithm https://en.wikipedia.org/wiki/Freivalds%27_algorithm#The_algorithm
- [9] Justin Thaler. [Proofs, Arguments, and Zero-Knowledge](https://people.cs.georgetown.edu/jthaler/ProofsArgsAndZK.pdf)