前言
众所周知,Python的for循环本质上要比C慢很多。 而且深度学习和机器学习算法严重依赖通过for循环执行的矩阵运算。
这就是为什么像numpy等这样包诞生,它们在numpy数组上提供向量化的操作。这意味着它将通常在Python中完成的for循环推进到C的级别。
我们希望将最大期望算法(Expectation-Maximization algorithm, EM)用于无监督学习(例如,识别MNIST数据集中的手写数字),并且我们的数据是二进制的(例如,二进制图像)。一种常见的方法是将数据建模为伯努利混合模型;一个人伯努利分布的加权和,如果每个分布有自己的标量权重π和自己的平均向量μ,并表示一组数据(例如,如果我们的数据是数字2、3&4的图形,我们使用3伯努利模型,一个伯努利将是数字2,另一个是4,等等)。总的来说,前者是一个向量,后者是一个矩阵。
Bernoulli mixture model
Distribution of one observation x given the cluster k
在E-step中,我们特别感兴趣的是隐变量后验的期望值,也就是所谓的responsibilities。
E-step of EM algorithm
γ实际返回的期望值观察n属于集群k。
γ是一个NxK矩阵;对于每个观测,我们分配的一个概率属于每个集群。最大值是我们指定的值。
因为:向量化过程中最重要的事情是要理解变量的维数。
Pipeline
我们将创建一个E_step函数来计算上面的表达式并用下面的代码进行测试:
第一次尝试
通过观察这些方程,我们可以看到有3个循环,每个例子 D 有一个循环,每个集群 K 有一个循环,每个对象 D 有一个循环,我们将按这个顺序循环。所以我们要每次用一个元素填充矩阵γ。
下图是结果:
第二次尝试
我们想去掉for loop D。因此,每个依赖于 D 的term应该变成一个向量。在for loop中,我们有两个变量;μ和x。因此 x 和 μ → 向量。问题是,它是 μ**x。
有一个函数,它把一个幂运算变成了乘法运算。没错,就是对数!因此,让我们使用对数来表示我们的表达式,然后对结果取指数。
关于对数概率的操作是首选的,因为它们提供了数值稳定性!
即使在我们的例子中它没有任何影响,每次你使用对数的时候,在表达式中使用一个常量 epsilon 来表示稳定性(不趋于0,是-inf)。
因此,我们将不得不对元素进行矢量乘法,easy!
结果是:
第三次尝试
一次一个loop:K turn
在向量化过程中,有如下操作:
标量→向量→矩阵
当我们用numpy数组替换越来越多的循环时,越来越多的代码将在C上运行。
我们使用之前的实现,我们想要删除K for loop。因此,每一个依赖于K的标量都会变成一个向量,每一个向量都会变成一个矩阵。这意味着X和μ将保持不变,π变成矩阵,γ变成向量。
结果:
n=1000的时候,我们只花了一半的时间!
第四次尝试
还有一个循环。我们可以有一个loop-python-free吗?come on!
由于我们要将矩阵*向量运算转换成矩阵@矩阵运算,我们需要取前者的传输矩阵(@是正则的矩阵乘法)。记住,现在我们的输出必须是整个γ矩阵。
一个循环也没有!代码看起来很优雅,只有三行!
对于n=1000,我们的运行时长从11.688下降到0.012!
总结
那么,当你想向量化一个表达式时,你需要做什么呢?