【程式碼閱讀】最大均值差異(Maximum Mean Discrepancy, MMD)損失函式程式碼解讀(Pytroch版)

【程式碼閱讀】最大均值差異(Maximum Mean Discrepancy, MMD)損失函式程式碼解讀(Pytroch版)
1 Star2 Stars3 Stars4 Stars5 Stars 給文章打分!
Loading...

程式碼及參考資料來源

Source code: easezyc/deep-transfer-learning [Github]
參考資料:遷移學習簡明手冊

MMD介紹

MMD(最大均值差異)是遷移學習,尤其是Domain adaptation (域適應)中使用最廣泛(目前)的一種損失函式,主要用來度量兩個不同但相關的分佈的距離。兩個分佈的距離定義為:

(1)MMD(X,Y)=||1n∑i=1nϕ(xi)−1m∑j=1mϕ(yj)||H2″ role=”presentation” style=”position: relative;”>MMD(X,Y)=||1n∑i=1nϕ(xi)−1m∑j=1mϕ(yj)||2H(1)(1)MMD(X,Y)=||1n∑i=1nϕ(xi)−1m∑j=1mϕ(yj)||H2

MMD(X,Y) = ||\frac{1}{n}\sum_{i=1}^n\phi(x_i)-\frac{1}{m}\sum_{j=1}^m\phi(y_j)||_H^2\tag{1}
其中 H” role=”presentation” style=”position: relative;”>HHH 表示這個距離是由 ϕ()” role=”presentation” style=”position: relative;”>ϕ()ϕ()\phi() 將資料對映到再生希爾伯特空間(RKHS)中進行度量的。

為什麼要用MMD?

Domain adaptation的目的是將源域(Source domain)中學到的知識可以應用到不同但相關的目標域(Target domain)。本質上是要找到一個變換函式,使得變換後的源域資料和目標域資料的距離是最小的。所以這其中就要涉及如何度量兩個域中資料分佈差異的問題,因此也就用到了MMD。至於Domain adaptation的前生今世可以參考王晉東大佬的知乎專欄

MMD的理論推導

MMD的關鍵在於如何找到一個合適的 ϕ()” role=”presentation” style=”position: relative;”>ϕ()ϕ()\phi() 來作為一個對映函式。但是這個對映函式可能在不同的任務中都不是固定的,並且這個對映可能高維空間中的對映,所以是很難去選取或者定義的。那如果不能知道ϕ” role=”presentation” style=”position: relative;”>ϕϕ\phi,那MMD該如何求呢?我們先展開把MMD展開:

(2)MMD(X,Y)=||1n2∑i,i′nϕ(xi)ϕ(xi′)−2nm∑i,jnϕ(xi)ϕ(yj) 1m2∑j,j′nϕ(yj)ϕ(yj′)||H2″ role=”presentation” style=”position: relative;”>MMD(X,Y)=||1n2∑i,i′nϕ(xi)ϕ(x′i)−2nm∑i,jnϕ(xi)ϕ(yj) 1m2∑j,j′nϕ(yj)ϕ(y′j)||2H(2)(2)MMD(X,Y)=||1n2∑i,i′nϕ(xi)ϕ(xi′)−2nm∑i,jnϕ(xi)ϕ(yj) 1m2∑j,j′nϕ(yj)ϕ(yj′)||H2

MMD(X,Y) =||\frac{1}{n^2}\sum_{i,i’}^n\phi(x_i)\phi(x_i’)-\frac{2}{nm}\sum_{i,j}^n\phi(x_i)\phi(y_j) \frac{1}{m^2}\sum_{j,j’}^n\phi(y_j)\phi(y_j’)||_H^2\tag{2}
展開後就出現了ϕ(xi)ϕ(xi′)” role=”presentation” style=”position: relative;”>ϕ(xi)ϕ(x′i)ϕ(xi)ϕ(xi′)\phi(x_i)\phi(x_i’)的形式,這樣聯絡SVM中的核函式k(∗)” role=”presentation” style=”position: relative;”>k(∗)k(∗)k(*),就可以跳過計算ϕ” role=”presentation” style=”position: relative;”>ϕϕ\phi的部分,直接求k(xi)k(xi′)” role=”presentation” style=”position: relative;”>k(xi)k(x′i)k(xi)k(xi′)k(x_i)k(x_i’)。所以MMD又可以表示為:

(3)MMD(X,Y)=||1n2∑i,i′nk(xi)k(xi′)−2nm∑i,jnk(xi)k(yj) 1m2∑j,j′nk(yj)k(yj′)||H2″ role=”presentation” style=”position: relative;”>MMD(X,Y)=||1n2∑i,i′nk(xi)k(x′i)−2nm∑i,jnk(xi)k(yj) 1m2∑j,j′nk(yj)k(y′j)||2H(3)(3)MMD(X,Y)=||1n2∑i,i′nk(xi)k(xi′)−2nm∑i,jnk(xi)k(yj) 1m2∑j,j′nk(yj)k(yj′)||H2

MMD(X,Y) =||\frac{1}{n^2}\sum_{i,i’}^nk(x_i)k(x_i’)-\frac{2}{nm}\sum_{i,j}^nk(x_i)k(y_j) \frac{1}{m^2}\sum_{j,j’}^nk(y_j)k(y_j’)||_H^2\tag{3}
在大多數論文中(比如DDC, DAN),都是用高斯核函式k(u,v)=e−||u−v||2σ” role=”presentation” style=”position: relative;”>k(u,v)=e−||u−v||2σk(u,v)=e−||u−v||2σk(u,v) = e^{\frac{-||u-v||^2}{\sigma}}來作為核函式,至於為什麼選用高斯核,最主要的應該是高斯核可以對映無窮維空間(具體的之後再分析)

理論到這裡就差不多了,那如何進行實現呢?

在TCA中,引入了一個核矩陣方便計算

(4)[Ks,sKs,sKs,tKt,t]” role=”presentation” style=”position: relative;”>[Ks,sKs,tKs,sKt,t](4)(4)[Ks,sKs,sKs,tKt,t]

\begin{bmatrix}
K_{s,s} & K_{s,s} \\
K_{s,t} & K_{t,t} \\
\end{bmatrix} \tag{4}

以及L矩陣:

(5)li,j={1/n2,xi,xj∈Ds1/m2,xi,xj∈Ds−1/nm,otherwise” role=”presentation” style=”position: relative;”>li,j=⎧⎩⎨⎪⎪1/n2,1/m2,−1/nm,xi,xj∈Dsxi,xj∈Dsotherwise(5)(5)li,j={1/n2,xi,xj∈Ds1/m2,xi,xj∈Ds−1/nm,otherwise

l_{i,j} =
\begin{cases}
1/{n^2}, & \text{$x_i, x_j\in D_s$} \\
1/{m^2}, & \text{$x_i, x_j\in D_s$} \\
-1/{nm},& \text{otherwise}
\end{cases} \tag{5}

在實際應用中,高斯核的σ” role=”presentation” style=”position: relative;”>σσ\sigma會取多個值,分別求核函式然後取和,作為最後的核函式。

程式碼解讀

import torch
def guassian_kernel(source, target, kernel_mul=2.0, kernel_num=5, fix_sigma=None):
'''
將源域資料和目標域資料轉化為核矩陣,即上文中的K
Params: 
source: 源域資料(n * len(x))
target: 目標域資料(m * len(y))
kernel_mul: 
kernel_num: 取不同高斯核的數量
fix_sigma: 不同高斯核的sigma值
Return:
sum(kernel_val): 多個核矩陣之和
'''
n_samples = int(source.size()[0]) int(target.size()[0])# 求矩陣的行數,一般source和target的尺度是一樣的,這樣便於計算
total = torch.cat([source, target], dim=0)#將source,target按列方向合併
#將total複製(n m)份
total0 = total.unsqueeze(0).expand(int(total.size(0)), int(total.size(0)), int(total.size(1)))
#將total的每一行都複製成(n m)行,即每個資料都擴充套件成(n m)份
total1 = total.unsqueeze(1).expand(int(total.size(0)), int(total.size(0)), int(total.size(1)))
#求任意兩個資料之間的和,得到的矩陣中座標(i,j)代表total中第i行資料和第j行資料之間的l2 distance(i==j時為0)
L2_distance = ((total0-total1)**2).sum(2) 
#調整高斯核函式的sigma值
if fix_sigma:
bandwidth = fix_sigma
else:
bandwidth = torch.sum(L2_distance.data) / (n_samples**2-n_samples)
#以fix_sigma為中值,以kernel_mul為倍數取kernel_num個bandwidth值(比如fix_sigma為1時,得到[0.25,0.5,1,2,4]
bandwidth /= kernel_mul ** (kernel_num // 2)
bandwidth_list = [bandwidth * (kernel_mul**i) for i in range(kernel_num)]
#高斯核函式的數學表示式
kernel_val = [torch.exp(-L2_distance / bandwidth_temp) for bandwidth_temp in bandwidth_list]
#得到最終的核矩陣
return sum(kernel_val)#/len(kernel_val)
def mmd_rbf(source, target, kernel_mul=2.0, kernel_num=5, fix_sigma=None):
'''
計算源域資料和目標域資料的MMD距離
Params: 
source: 源域資料(n * len(x))
target: 目標域資料(m * len(y))
kernel_mul: 
kernel_num: 取不同高斯核的數量
fix_sigma: 不同高斯核的sigma值
Return:
loss: MMD loss
'''
batch_size = int(source.size()[0])#一般預設為源域和目標域的batchsize相同
kernels = guassian_kernel(source, target,
kernel_mul=kernel_mul, kernel_num=kernel_num, fix_sigma=fix_sigma)
#根據式(3)將核矩陣分成4部分
XX = kernels[:batch_size, :batch_size]
YY = kernels[batch_size:, batch_size:]
XY = kernels[:batch_size, batch_size:]
YX = kernels[batch_size:, :batch_size]
loss = torch.mean(XX   YY - XY -YX)
return loss#因為一般都是n==m,所以L矩陣一般不加入計算

程式碼示例

為了體現以上程式碼的有效性,我們參考連結生成了兩組不同分佈的資料。

import random
import matplotlib
import matplotlib.pyplot as plt
SAMPLE_SIZE = 500
buckets = 50
#第一種分佈:對數正態分佈,得到一箇中值為mu,標準差為sigma的正態分佈。mu可以取任何值,sigma必須大於零。
plt.subplot(1,2,1)
plt.xlabel("random.lognormalvariate")
mu = -0.6
sigma = 0.15#將輸出資料限制到0-1之間
res1 = [random.lognormvariate(mu, sigma) for _ in xrange(1, SAMPLE_SIZE)]
plt.hist(res1, buckets)
#第二種分佈:beta分佈。引數的條件是alpha 和 beta 都要大於0, 返回值在0~1之間。
plt.subplot(1,2,2)
plt.xlabel("random.betavariate")
alpha = 1
beta = 10
res2 = [random.betavariate(alpha, beta) for _ in xrange(1, SAMPLE_SIZE)]
plt.hist(res2, buckets)
plt.savefig('data.jpg)
plt.show()

兩種資料分佈如下圖
這裡寫圖片描述
兩種分佈有明顯的差異,下面從兩個方面用MMD來量化這種差異:
1. 分別從不同分佈取兩組資料(每組為10*500)

from torch.autograd import Variable
#引數值見上段程式碼
#分別從對數正態分佈和beta分佈取兩組資料
diff_1 = []
for i in range(10):
diff_1.append([random.lognormvariate(mu, sigma) for _ in xrange(1, SAMPLE_SIZE)])
diff_2 = []
for i in range(10):
diff_2.append([random.betavariate(alpha, beta) for _ in xrange(1, SAMPLE_SIZE)])
X = torch.Tensor(diff_1)
Y = torch.Tensor(diff_2)
X,Y = Variable(X), Variable(Y)
print mmd_rbf(X,Y)

輸出結果為

Variable containing:
6.1926
[torch.FloatTensor of size 1]

2. 分別從相同分佈取兩組資料(每組為10*500)

from torch.autograd import Variable
#引數值見以上程式碼
#從對數正態分佈取兩組資料
same_1 = []
for i in range(10):
same_1.append([random.lognormvariate(mu, sigma) for _ in xrange(1, SAMPLE_SIZE)])
same_2 = []
for i in range(10):
same_2.append([random.lognormvariate(mu, sigma) for _ in xrange(1, SAMPLE_SIZE)])
X = torch.Tensor(same_1)
Y = torch.Tensor(same_2)
X,Y = Variable(X), Variable(Y)
print mmd_rbf(X,Y)

輸出結果為

Variable containing:
0.6014
[torch.FloatTensor of size 1]

可以明顯看出同分布資料和不同分佈資料之間的差距被量化了出來,且符合之前理論所說:不同分佈MMD的值大於相同分佈MMD的值。
PS,在實驗中發現一個問題,就是取資料時要在0-1的範圍內取,不然MMD就失效了。
如果錯誤之處,請指正,感謝閱讀

相關文章

程式語言 最新文章