完全解析：使用Faiss进行海量特征的相似度匹配

背景

我们不妨想象下面的几个例子：

1. 输入一张商品的图片，从商品库中匹配出相似的商品，这是以图搜图的一个例子；
2. 输入一小段音乐，从音乐库中匹配出对应的音乐出，这是MIR的一个例子；
3. 输入一张人脸，从人脸底库中匹配出对应的人，这是1:N 人脸识别的一个例子；

像这样的例子还有很多，事实上，以神经网络对样本进行特征的提取，然后在海量的特征库里进行特征相似度的搜索/比对/匹配，已经是AI技术落地的一大领域。Faiss就是Facebook维护的一个高效的特征相似度匹配和聚类的库。

本文将从最基本的特征比对说起，然后落脚到我们为什么需要Faiss，以及Faiss上提供的在特征比对之外的功能。

一个简单的特征相似度比对的例子

设想我们使用一个在ImageNet上预训练的resnet50模型来提特征，因为只需要最后的2048维特征，我们在例子中把resnet50网络最后的fc层去掉。别着急，在deepvac项目的examples目录下正好有一个完全一样的例子：

https://github.com/DeepVAC/deepvac/blob/master/examples/a_resnet_project/test_emb.py

假设我们现在要在db里放入7030张图片的特征来作为我们的特征库，之后，待搜索的图片就和该特征库来做相似度匹配。为了实验这个想法，首先我们来制作特征库。当我们执行程序test_emb.py后，在process函数中，该test_emb.py会将目录下的7030张图片提取为7030个2048维的特征：

2020-09-02 11:50:44,507 1073:master INFO     feature shape: torch.Size([1, 2048])
2020-09-02 11:50:44,508 1073:master INFO     feature db shape: torch.Size([1, 2048])
2020-09-02 11:50:44,527 1073:master INFO     feature shape: torch.Size([1, 2048])
2020-09-02 11:50:44,529 1073:master INFO     feature db shape: torch.Size([2, 2048])
2020-09-02 11:50:44,544 1073:master INFO     feature shape: torch.Size([1, 2048])
2020-09-02 11:50:44,544 1073:master INFO     feature db shape: torch.Size([3, 2048])
2020-09-02 11:50:44,571 1073:master INFO     feature shape: torch.Size([1, 2048])
2020-09-02 11:50:44,571 1073:master INFO     feature db shape: torch.Size([4, 2048])
2020-09-02 11:50:44,599 1073:master INFO     feature shape: torch.Size([1, 2048])
2020-09-02 11:50:44,600 1073:master INFO     feature db shape: torch.Size([5, 2048])
......

那么最后这个特征库的shape就为7030 x 2048。现在我们要检索一张图片该怎么做呢？对该图片提特征、和db里的特征进行L2距离的计算、找出距离最近的1个或k个。代码如下所示：

#deepvac关于配置文件的标准
>>> from config import config

#需要使用test_emb.py中的DeepvacEmb类
>>> from test_emb import DeepvacEmb
>>> emb = DeepvacEmb(config)

#加载上文保存的特征库
>>> emb.loadDB('resnet50_gemfield_test.emb')

#获得测试图片的2048维特征
>>> xq = emb.getPillowImgFeatureVector("/gemfield/hostpv/nsfw/test/porn/00002640.000000.jpg")
>>> xq.shape
torch.Size([1, 2048])

#从特征库中搜索最匹配的，D为距离，I为索引
>>> D,I = emb.search(xq)
>>> D
[0.0]
>>> I
[248]

#从特征库中搜索最匹配的3个
>>> D,I = emb.search(xq,k=3)
>>> D
[0.0, 14.658945083618164, 15.17888069152832]
>>> I
[248, 225, 448]

上面的代码演示了从特征库中去搜索最相似的图片。其中使用到的Deepvac的search API就是基于PyTorch的torch.norm() API进行的L2距离的计算。

在多个CUDA设备上进行特征的相似度搜索

有的时候，我们的特征库太大，以至于一个CUDA卡已经装载不下（比如20G的特征库，而一个RTX2080ti才11G显存）。这个时候该怎么办呢？一个直觉的想法就是将特征库拆分成多份，分别放在多个CUDA设备上。然后我们的查询向量xq在每个CUDA设备上和特征库相比对，最后再汇总排序取出距离最小的那个。

deepvac的syszux_feature_vector模块就是这么干的，你可以看下该模块中的NamesPathsClsFeatureVector类的实现：

https://github.com/DeepVAC/deepvac/blob/master/deepvac/syszux_feature_vector.pygithub.com

这个类也是DeepVAC内部常用的一个工具类。现在问题来了，即使有Deepvac中的search()方法，即使还有syszux_feature_vector模块中的searchAllDB()方法，但在以下方面还是有点不够灵活和有力：

1. 如何一次搜索一批的特征，而不仅仅是一个特征？
2. 如何返回更相似度最近的一批特征，而不只是一个特征？（好吧，Deepvac类也支持）
3. 如何让特征库使用的内存空间更小？（你看，上面都需要把特征库拆分到多个cuda设备上了）
4. 搜索速度方面如何更快？

这就是Faiss库存在的意义。Faiss：Facebook AI Similarity Search。

Faiss环境准备

提前说明的福利：你可以使用如下的docker环境，从而省却自己配置环境的烦恼：

#只使用cpu
docker run -it gemfield/faiss:1.6.3-devel bash

#使用GPU的话
docker run --gpus all -it gemfield/faiss:1.6.3-devel bash

如果你不想使用上述的Docker镜像，那么需要自行安装依赖、编译、安装，这些步骤至少有：

1，安装依赖包

安装libopenblas-dev：

apt install libopenblas-dev

否则报错：CMake Error at /usr/share/cmake-3.18/Modules/FindPackageHandleStandardArgs.cmake:165 (message):
Could NOT find BLAS (missing: BLAS_LIBRARIES)

安装swig：

apt install swig

否则报错：Could NOT find SWIG (missing: SWIG_EXECUTABLE SWIG_DIR python)。

2，编译

在Faiss目录里，执行如下操作：

make build
cmake ..
make VERBOSE=1
make install

整个编译的最终产物有：

• libfaiss.a （C++库，也用于链接_swigfaiss.so）
• libfaiss_python_callbacks.a （用于链接_swigfaiss.so）
• _swigfaiss.so（用于python）

最后一步的make install会将libfaiss.a安装到/usr/local/lib/libfaiss.a，以及将头文件安装到/usr/local/include/faiss/目录下。

3，安装python

在build/faiss/python目录下：

python setup.py install

会将把_swigfaiss.so安装在/opt/conda/lib/python3.7/site-packages/faiss-1.6.3-py3.7.egg/faiss/

Faiss的简单使用：Flat

我们先定义两个变量xb和xq。xb用来表示特征库，xq用来表示待查询的2048维向量。如果沿用上面的例子，则xb就是提前存储了7030个样本的特征的“数据库”，它的shape就是7030x2048——这样的“数据库”在Faiss中称作Index object。

就像数据库的种类有很多一样，Index的种类也有很多，我们先从最简单的IndexFlatL2开始。IndexFlatL2是什么类型的“数据库”呢？就是使用暴力L2搜索的数据库——也就是和特征库中的每个特征进行L2距离计算然后取出距离最近的那个。是不是看着很熟悉？没错，这和上文中提到的DeepVAC的search() API的原理是一模一样的。

好多种Index在被真正检索前，是需要一个“训练”操作的，类似给数据库的字段建立索引，Index object的“训练”其实就是提前根据特征的分布进行聚类训练。而我们的IndexFlatL2并不在列，前面说了，它是最简单的“数据库”：

import numpy as np
import faiss

d = 2048                          # dimension
nb = 7030                     # database size
nq = 10                       # nb of queries

np.random.seed(1234)             # make reproducible
xb = np.random.random((nb, d)).astype('float32')
xb[:, 0] += np.arange(nb) / 1000.

xq = np.random.random((nq, d)).astype('float32')
xq[:, 0] += np.arange(nq) / 1000.

index = faiss.IndexFlatL2(d)
print(index.is_trained)

index.add(xb)
print(index.ntotal)

k = 4                          # we want to see 4 nearest neighbors
D, I = index.search(xb[:5], k) # sanity check
print("I: ",I)
print("D: ",D)

上面是使用IndexFlatL2的一个简单例子，Gemfield初始化了一个7030x2048的特征库，然后一次检索5个特征，也就是xq是5x2048。代码输出：

True
7030
I:  [[   0  998 1309 1134]
 [   1  252 2013 2790]
 [   2 1046 1342 1204]
 [   3   56 1846  910]
 [   4  714  324   50]]
D:  [[  0.      323.07217 325.1012  325.8144 ]
 [  0.      308.83826 313.52512 317.5051 ]
 [  0.      313.24506 313.52145 316.1954 ]
 [  0.      311.7074  313.01157 313.31152]
 [  0.      316.34128 316.5642  317.65128]]

I是5个待检索特征各自前4个特征库中最相似的index，D是对应的距离。既然前文说过，IndexFlatL2只是最简单的一种“数据库”，那么其它复杂的数据库都有什么呢？你可以使用如下的方法简单list下就可以管中窥豹了：

>>> import faiss
>>> for s in dir(faiss):
...   if s.startswith('Index'):
...     print(s)
... 
IndexBinary
IndexBinaryFlat
IndexFlat
IndexFlat1D
IndexFlat1D_swigregister
IndexFlatIP
IndexFlatIP_swigregister
IndexFlatL2
IndexFlatL2BaseShift
IndexFlatL2BaseShift_swigregister
IndexFlatL2_swigregister
IndexFlat_swigregister
IndexHNSW
IndexHNSW2Level
IndexHNSW2Level_swigregister
IndexHNSWFlat
IndexHNSWFlat_swigregister
IndexHNSWPQ
IndexHNSWPQ_swigregister
IndexHNSWSQ
IndexHNSWSQ_swigregister
IndexHNSW_swigregister
IndexIDMap
IndexIDMap2
IndexIDMap2_swigregister
IndexIDMap_swigregister
IndexIVF
IndexIVFFlat
IndexIVFFlatDedup
IndexIVFFlatDedup_swigregister
IndexIVFFlat_swigregister
IndexIVFPQ
IndexIVFPQR
IndexIVFPQR_swigregister
IndexIVFPQStats
IndexIVFPQStats_swigregister
IndexIVFPQ_swigregister
IndexIVFScalarQuantizer
IndexIVFScalarQuantizer_swigregister
IndexIVFSpectralHash
IndexIVFSpectralHash_swigregister
IndexIVFStats
IndexIVFStats_swigregister
IndexIVF_swigregister
IndexLSH
IndexLSH_swigregister
IndexLattice
IndexLattice_swigregister
IndexPQ
IndexPQStats
IndexPQStats_swigregister
IndexPQ_swigregister
......

以及GPU上的Index类型：

>>> import faiss
>>> for s in dir(faiss):
...   if s.startswith('GpuIndex'):
...     print(s)
... 
GpuIndexFlat
GpuIndexFlatConfig
GpuIndexFlatConfig_swigregister
GpuIndexFlatIP
GpuIndexFlatIP_swigregister
GpuIndexFlatL2
GpuIndexFlatL2_swigregister
GpuIndexFlat_swigregister
GpuIndexIVF
GpuIndexIVFConfig
GpuIndexIVFConfig_swigregister
GpuIndexIVFFlat
GpuIndexIVFFlatConfig
GpuIndexIVFFlatConfig_swigregister
GpuIndexIVFFlat_swigregister
GpuIndexIVFPQ
GpuIndexIVFPQConfig
GpuIndexIVFPQConfig_swigregister
GpuIndexIVFPQ_swigregister
GpuIndexIVFScalarQuantizer
GpuIndexIVFScalarQuantizerConfig
GpuIndexIVFScalarQuantizerConfig_swigregister
GpuIndexIVFScalarQuantizer_swigregister
GpuIndexIVF_swigregister
GpuIndex_swigregister

这里已经看到了各种CPU的Index和GPU的Index，那么实践中这两者谁更快呢？

• xq的batch很小，Index很小：CPU通常更快；
• xq的batch很小，Index很大：GPU通常更快；
• xq的batch很大，Index很小：随便；
• xq的batch很大，Index很大：GPU通常更快；
• GPU通常比CPU快5到10倍；

让Faiss使用更少的内存：PQ

IndexFlatL2的暴力L2距离匹配是最基本的用法。很快你就会遇到一个问题，当特征库很大的时候，一个显存就很快装载不下特征库了；即使是换成使用内存来装载特征库，也装载不下了。有没有一种压缩算法能减小特征库的大小呢？这就是PQ（Product Quantizer）。

不管是IndexFlatL2、IndexIVFFlat、或者名字中包含有Flat的Index，都会在加载到内存中的特征库中保存全量的特征，以2048维的float类型为例，一个样本就需要8192个字节。如果类似天猫京东这样有1亿件商品，每件商品就按照一个样本图片来算，也需要...763G的特征库。如果我们能将特征库进行有效压缩，那可以节省相当可观的资源。PQ就是带着这样的目的而来的，这是一种有损压缩，所以这种Index进行检索的返回值只是近似的准确，而不是像IndexFlatL2那样可以返回绝对准确的值。那么PQ是如何将一个样本的8192字节压缩的更小的呢？假设我们有100万个样本，每个样本有2048维特征（100万 x 2048）：

1. 我们将每个样本的2048维vector拆分为8个sub-vector，每个sub-vector就是256维了。这样就会有8个100万x256维的矩阵；
2. 我们在这8个矩阵上使用k = 256的k-means 聚类算法（Gemfield：这里的256和上面的256没啥关系），这样每个矩阵上会得到256个centroid（质心），一共有8组256个centroid；聚类算法的使用，使得每个256维的centroid成为最能代表那3906个256维特征的一个向量（100w/256 = 3906）；为啥选择k=256呢？因为256个centroid的索引可以使用一个字节装下；
3. 我们再回到那1亿个商品的特征库，我们将每个商品的2048维特征拆分成8个sub-vector，每一个sub-vector都被替换为距离最近的那个centroid，并使用该centroid的id来表示（1个字节！）。如此以来，2048维特征就被替换成了8个centroid ID，也就是8个字节！如此以来，商品特征库就从763G下降到了0.745G！

内存的使用量确实降下来了，但是如果特征库只包含centroid ID的话，怎么进行向量的相似度计算呢？只有centroid ID的话，怎么计算L2距离呢？？？是不是需要把特征库里的1亿个特征（每个特征此刻为8个字节）中的centroid ID替换成centroid来进行L2距离的计算呢？这固然是一个直截了当的想法，只不过，替换回去（也就是解压缩）的话，每个特征本来用8个字节表示，如今又要变成8096个字节了，那我们这么折腾一回图啥呀！相反，我们可以这么做：

1，一个2048维的xq来了，将该xq vector拆分为8个sub-vector；

2，xq vector的每个sub-vector就和8x256的centroid中对应的那组1x256的centroid进行L2距离计算，得到256个距离；那么8组下来，就会得到8组256个L2距离；这个计算量仿佛和一个只有256个样本特征的特征库差不多；

3，上述的8组256个L2的距离，可以看成是一个表：256行x8列（毕竟是一个样本的特征被拆分了8个sub vector，所以不是8x256而是256x8）；这个表很重要啊，我暂且称之为gemfield表；

4，特征库现在是已经压缩的状态，也就是说每个样本特征用的是8个centroid ID来表示的。而每个特征库的样本（每个xb记录）和xq的L2距离，就等于每个xb记录的8个sub-vector和xq的8个sub-vector的L2距离之和，就约等于xb记录的8个centroid（中的每个）和xq的8个sub-vector（中的每个）的L2距离之和！而xb记录的8个centorid中的每个和xq的8个sub-vector中的每个之间的L2距离，直接通过上述的gemfield表查表获得！

5，将距离排序，取得距离最近的k个xb记录（的index和distance）。

上面的centroid（最类中心的那个2048维向量），在Faiss中我们称之为code；上述的8组256个centroid，在Faiss中我们称之为8个code book；这8个code book可以表达256^8个值，还是很大的。

让Faiss进行更快的检索：IVF

IndexFlatL2的暴力L2距离匹配是最基本的用法。很快你又会遇到一个问题，当检索量很大的时候，每次检索哪怕减少一点时间，对整体系统的性能提升也会有很大的帮助；换言之，可以提升边际效益。那么如何让检索更快呢？事实上，更快的检索来自于两个方面：

1. 两两特征比对更少的计算量；PQ顺带着做了；
2. 只和特征库的一部分进行比对；和特征库的每一个特征进行比对，叫做穷举；只和部分特征进行比对，叫做IVF；

问题是，为什么和特征库的一部分进行比对就能找到想要的答案呢？呃，不能。只能找到近似正确的答案。为什么和特征库的一部分进行比对就能找到近似正确的答案呢？呃，倒排索引（IVF）。IVF用来提前过滤dataset从而避开对全部特征库的穷举搜索，比如gemfield使用k-means算法将数据库进行聚类（比如100个分区），当查询一个xq的时候，和能代表这100个分区的100个centroid进行L2比较，拿到最接近的5个；然后在这5个分区里再进行穷举搜索。这样，有95%的xb记录就被忽略掉了。

IVF，全称Inverted File Index，中文翻译为倒排索引，指的是在文本搜索中（比如搜索引擎），将每个单词到该单词所属的文档之间的映射关系保存到数据库中。Gemfield是非常不喜欢这个概念的，英文的Inverted File Index我都感觉词不达意，而中文的“倒排索引”更是让人摸不着头脑。以一本书为例，书前面的目录就是索引：

index，索引

书最后面的index就是倒排索引，如下所示：

IVF，倒排索引

因此，在Faiss中所有带有IVF的index，指的都是提前将数据库使用k-means聚类算法划分为多个partition，每个partition中包含有对应的feature vectors（这就是inverted file lists指向的），同时，每个partition还有对应的centroid，用来决定该partition是应该被穷举搜索，还是被忽略。

IVF这个方法当然加快了检索的速度，但是也是有代价的。假如xq是落在partition的边缘地带，那这个xq最近的记录大概率是距离它最近（这个最近指的是和centroid比较）的前面好多个partition里的一个，这样即使我们指定很多个比较近的partition（比如5个），如果ground truth实际上是在第6个partition里，这样返回的结果依然是不准确的。

因此，带有IVF的检索只能返回近似正确的值，而不能像IndexFlatL2那样返回绝对正确的值。但是像上面这种ground truth是在第6个partition里，而我们指定5个partition也是比指定1个partition要更准确些的（代价就是检索时间多了）。一旦找到最近的那个partition后，便就要开始穷举检索该partition里的每条记录了，这就是IndexIVFFlat的由来。在某个partition中进行搜索的过程还可以使用上一节的PQ压缩的算法，因此，在Faiss中，我们还经常会使用的一个Index叫作IndexIVFPQ。

这个partition的概念，在Faiss中称之为Voronoi cells；选择某个Voronoi cells然后进行检索的动作，称之为“probe”；而在最近的“多少个”Voronoi cells中进行检索，这个“多少个”的数量称之为nprobe；在IVF Index object的属性中，就有nprobe这个属性：

import numpy as np
import faiss

d = 2048                          # dimension
nb = 7030                     # database size
nq = 10                       # nb of queries

np.random.seed(1234)             # make reproducible
xb = np.random.random((nb, d)).astype('float32')
xb[:, 0] += np.arange(nb) / 1000.

xq = np.random.random((nq, d)).astype('float32')
xq[:, 0] += np.arange(nq) / 1000.

index = faiss.IndexFlatL2(d)
print(index.is_trained)

index.add(xb)
print(index.ntotal)

######
nlist = 100
m = 8
k = 4
quantizer = faiss.IndexFlatL2(d)  # this remains the same
index = faiss.IndexIVFPQ(quantizer, d, nlist, m, 8)
                                  # 8 specifies that each sub-vector is encoded as 8 bits
index.train(xb)
index.add(xb)
#nprobe
index.nprobe = 10
D, I = index.search(xb[:5], k)
print("I: ",I)

IndexIVFPQ的参数中，d代表特征的维度2048，nlist代表Voronoi cells的数量，m代表subquantizers的数量（给PQ用），8代表使用8bits表示一个特征的1个sub-vector。代码输出：

root@pytorch160-ai2:/home/gemfield# python test.py 
#当index.nprobe = 10
I:  [[  14   68    0  599]
 [ 297    1  546  227]
 [ 246  311    2  137]
 [ 213    3  256  248]
 [   4  919  644 1351]]

可见确实只能返回近似准确的值。

无论如何内存都不够用：Distributed index

在某些情况下，内存或者显存无论如何都不够使用。比如要加载20多个GB的特征库，但是显卡只有11GB，又不想妥协任何的准确度。那怎么办呢？方法之一就是将特征库分散加载到多个服务器的CUDA设备上。

在项目仓库的faiss/contrib/目录下，Faiss提供了rpc.py，一个小型的RPC库。基于该RPC库，Faiss在faiss/contrib/client_server.py中实现了run_index_server() API和ClientIndex类。在不同的服务器上，使用run_index_server() API运行服务：

#比如四个服务器: ai{1..4}.gemfield.org
#index就是所有特征库的四分之一
run_index_server(index, port, v6=v6)

在客户端，使用client_index=ClientIndex(host_port_array)来访问：

syszux_ports = [
    ('ai1.gemfield.org', 12010),
    ('ai2.gemfield.org', 12011),
    ('ai3.gemfield.org', 12012),
    ('ai4.gemfield.org', 12013),
]

client_index = ClientIndex(syszux_ports)
#xq是query vector
D, I = client_index.search(xq, 5)

其实client_index获得分段的结果，最后归并取最好结果的这个过程，和文章开头提到的syszux_feature_vector.py里的逻辑是一样的。

无论如何内存都不够用：On-disk storage

参考项目中的faiss/contrib/ondisk.py。

当xq是pytorch的tensor时

在前文的各种例子中，你都发现，无论是xb还是xq，它们都是转换为numpy数组才开始调用Faiss的API的。没错，在Faiss中，numpy就是一等公民。既然Faiss是Facebook的项目，既然Facebook还有PyTorch这么流行的项目，那么在Faiss中为什么不原生支持PyTorch的Tensor呢？更何况numpy只能存活在CPU上，Tensor可是CPU和GPU通吃啊！实在是令人莫名其妙、瞠目结舌、匪夷所思！

在前文的例子中：

https://github.com/DeepVAC/deepvac/blob/master/examples/a_resnet_project/test_emb.pygithub.com

我们的特征库可都是使用PyTorch的Tensor来存储和序列化的，查询特征xq也是tensor，总不能每次都从Tensor转换成numpy吧。目前，Faiss提供了一种临时的措施，可以直接读取Tensor底层的内存（Tensor必须是is_contiguous的），然后使用index的search_c() API来检索。关于在index中检索PyTorch Tensor的详细用法，你可以参考syszux_feature_vector.py中的NamesPathsClsFeatureVectorByFaissPytorch类：

https://github.com/DeepVAC/deepvac/blob/master/deepvac/syszux_feature_vector.pygithub.com

最后，如何选择一种Index呢？

我们已经见识过的关键字有Flat、IVF、PQ，那么如何选择一种Index来匹配我们的场景呢？

1. 当需要绝对准确的结果时，使用Flat；比如IndexFlatL2 或者 IndexFlatIP；
2. 如果内存完全够用富裕的不行，使用HNSW；如果一般够用，使用Flat；如果有点吃紧，使用PCARx,...,SQ8；如果非常吃紧，使用OPQx_y,...,PQx；
3. 如果特征库小于1M个记录，使用"...,IVFx,..."；如果在1M到10M之间，使用"...,IVF65536_HNSW32,..."；如果在10M - 100M之间，使用"...,IVF262144_HNSW32,..."；如果在100M - 1B之间，使用"...,IVF1048576_HNSW32,..."。

用于训练的xb记录越多，聚类算法的训练效果越好，但训练需要花的时间也就越多，别忘了这一点。