UKRの実装と理解@多賀

ukr/sample/ukr_gaussian.py

@@ -95,7 +95,7 @@ def estimate_e(self, X, Y, Z, R):
                                d_in,
                                δ_in,
                                optimize=True)
-        diff = 2 * (diff_left - diff_right) / X.shape[0]
+        diff = 2 * (diff_left - diff_right) / (X.shape[0] * self.σ**2)
         Z -= self.η * (diff + 2 * self.λ * Z / X.shape[0])
         return Z
 
@@ -105,7 +105,7 @@ def estimate_e(self, X, Y, Z, R):
     from visualizer import visualize_history
 
 
-    X = gen_saddle_shape(num_samples=200, random_seed=1, noise_scale=0.001)
-    ukr = UKR(latent_dim=2, eta=10, rambda=1e-1, sigma=1, scale=1e-2)
-    history = ukr.fit(X, num_epoch=200)
+    X = gen_saddle_shape(num_samples=400, random_seed=1, noise_scale=0.05)
+    ukr = UKR(latent_dim=1, eta=5, rambda=0, sigma=0.1, scale=1e-2)
+    history = ukr.fit(X, num_epoch=100)
     visualize_history(X, history['f'], history['Z'], save_gif=False)

ukr/taga/cross_test_taga.py

@@ -0,0 +1,39 @@
+import numpy as np
+import data
+from ukr_gaussian import UKR as GaussianUKR
+from ukr import UKR as tagaUKR
+from matplotlib import pyplot as plt
+
+seed = 1
+
+if __name__ == '__main__':
+    N = 400
+    X = data.gen_saddle_shape(num_samples=N, random_seed=seed, noise_scale=0.05)
+
+    params_for_gukr = dict(
+        latent_dim=2,
+        eta=1,
+        rambda=0,
+        sigma=0.1,
+        scale=1e-2
+    )
+    params_for_tagaukr = dict(
+        N=N,
+        D=3,
+        L=2,
+        eta=1,
+        sigma=0.1,
+        scale=1e-2
+    )
+
+    ukr_1 = GaussianUKR(**params_for_gukr)
+    history_1 = ukr_1.fit(X, num_epoch=100, seed=seed, f_resolution=10)
+    ukr_2 = tagaUKR(**params_for_tagaukr)
+    history_2 = ukr_2.fit(X, num_epoch=100, seed=seed, resolution=10)
+
+    is_Y_close = np.allclose(history_1['Y'], history_2['Y'])
+    is_f_close = np.allclose(history_1['f'], history_2['f'])
+    is_Z_close = np.allclose(history_1['Z'], history_2['Z'])
+    print("Y: ", is_Y_close)
+    print("f: ", is_f_close)
+    print("Z: ", is_Z_close)

ukr/taga/data.py

@@ -28,6 +28,22 @@ def gen_2d_sin_curve(num_samples, random_seed=None, noise_scale=0.01):
     return X
 
 
+def gen_helix(num_samples, random_seed=None, noise_scale=0.05):
+    np.random.seed(random_seed)
+    z = np.random.uniform(low=-1, high=+1, size=(num_samples,))
+
+    pi = np.pi
+    t = z * pi * 3
+
+    X = np.empty((num_samples,3))
+    X[:, 0] = np.cos(t)
+    X[:, 1] = np.sin(t)
+    X[:, 2] = t
+    X += np.random.normal(loc=0, scale=noise_scale, size=X.shape)
+
+    return X
+
+
 if __name__ == '__main__':
     from matplotlib import pyplot as plt
 

ukr/taga/ukr.py

@@ -1,36 +1,100 @@
 import numpy as np
 from tqdm import tqdm
+from utils import make_grid
 
+seed = 1
 
 class UKR(object):
-    def __init__(self, _):
-        NotImplemented
+    def __init__(self, N, D, L, eta, sigma, scale):
+        self.N = N
+        self.D = D
+        self.L = L
+        self.eta = eta
+        self.sigma = sigma
+        self.scale = scale
+
+    #データの距離
+    def distance(self, A, B):
+        #Aはデータ(こっちで用意するやつ),Bは構造(ノード)(NxN)
+        Dist = np.sum((A[:, None, :]-B[None, :, :])**2, axis=2)
+        return Dist
+
+    #写像の推定
+    def estimate_y(self, Z1, Z2, X):
+        #ガウス関数による重み行列k(z_n,z_i)(NxN)
+        k = np.exp((-1/(2*(self.sigma**2)))*self.distance(Z1, Z2))
+        #kをiでsumした行列K(z_n)(Nx1)
+        K = np.sum(k, axis=1) 
+        Y = (k@X)/K[:,None]
+        return Y
+
+    #潜在変数の推定
+    def estimate_z(self, Z1, Z2, X, Y):
+        #ガウス関数による重み行列k(z_n,z_i)(NxN)
+        k = np.exp((-1/(2*(self.sigma**2)))*self.distance(Z1, Z2))
+        #kをiでsumした行列K(z_n)(Nx1)
+        K = np.sum(k, axis=1)
+
+        r = k/K[:,None]
+        d = Y[:,None,:] - X[None,:,:]
+        delta = Z1[:,None,:] - Z2[None,:,:]
+
+        B = Y - X
+        A = np.sum(B[:,None,:]*d,axis=2)
+
+        C = r*A
+        CC = C+C.T
+        E_bibun = (2/(self.N*(self.sigma**2)))*np.sum(CC[:,:,None]*delta,axis=1)
+
+        Z_new = Z1 - self.eta*E_bibun
+
+        return Z_new
+
+    def fit(self, X, num_epoch, seed, resolution):
+        np.random.seed(seed)
+        Z = np.random.normal(scale=self.scale, size=(self.N, self.L))
 
-    def fit(self, _):
         history = dict(
-            E=np.zeros(()),
-            Y=np.zeros(()),
-            f=np.zeros(()),
-            Z=np.zeros(())
+            Y=np.zeros((num_epoch, self.N, self.D)),
+            f=np.zeros((num_epoch, resolution**self.L, self.D)),
+            Z=np.zeros((num_epoch, self.N, self.L))
         )
 
-        for epoch in tqdm(range(0)):
-            ...
+        for epoch in tqdm(range(num_epoch)):
+            Y = self.estimate_y(Z, Z, X)
+            Z = self.estimate_z(Z, Z, X, Y)
+
+            """
+            A = np.linspace(Z[:,0].min(), Z[:,0].max(), resolution)
+            B = np.linspace(Z[:,1].min(), Z[:,1].max(), resolution)
+            XX, YY = np.meshgrid(A,B)
+            M = np.concatenate([XX.reshape(-1)[:,None], YY.reshape(-1)[:,None]], axis=1)
+            """
+            #sampleプログラムに合わせた初期値，make_gridの中身はよくわかってない
+            Z_new = make_grid(resolution,
+                              bounds=(np.min(Z), np.max(Z)),
+                              dim=self.L)
+            f = self.estimate_y(Z_new, Z, X)
 
             history['Y'][epoch] = Y
             history['f'][epoch] = f
             history['Z'][epoch] = Z
-            history['E'][epoch] = ...
+
         return history
 
 
 if __name__ == '__main__':
     import data
     from visualizer import visualize_history
+    N = 400
+    D = 3
+    L = 2
+    T = 100
+    eta = 1
+    sigma = 0.1
+    resolution = 20
 
-    # 動かなさそうで動いてしまうプログラム
-    # お好きに書き直してください
-    X = data.gen_saddle_shape(num_samples=200, random_seed=0, noise_scale=0.05)
-    ukr = UKR(...)
-    history = ukr.fit(...)
-    visualize_history(X, history, save_gif=False)
+    X = data.gen_saddle_shape(num_samples=400, random_seed=seed, noise_scale=0.05)
+    ukr = UKR(N, D, L, eta, sigma, scale=1e-2)
+    history = ukr.fit(X, num_epoch=T, seed=seed, resolution=resolution)
+    visualize_history(X, history['Z'], history['f'], colormap=X[:,0], resolution=20, T=100, save_gif=False, filename="tmp")

ukr/taga/ukr_gauss.py

@@ -0,0 +1,114 @@
+from typing import OrderedDict
+import numpy as np
+from tqdm import tqdm
+from utils import make_grid
+from collections import OrderedDict
+
+seed = 1
+
+class UKR(object):
+    def __init__(self, N, D, L, eta, sigma, rambda, scale):
+        self.N = N
+        self.D = D
+        self.L = L
+        self.eta = eta
+        self.sigma = sigma
+        self.rambda = rambda
+        self.scale = scale
+
+    #データの距離
+    def distance(self, A, B):
+        #Aはデータ(こっちで用意するやつ),Bは構造(ノード)
+        Dist = np.sum((A[:, None]-B[None, :])**2, axis=2)
+        return Dist
+
+    #写像の推定
+    def estimate_y(self, Z1, Z2, X):
+        #ガウス関数による重み行列k(z_n,z_i)(NxN)
+        k = np.exp((-1/(2*(self.sigma**2)))*self.distance(Z1, Z2))
+        #kをiでsumした行列K(z_n)(Nx1)
+        K = np.sum(k, axis=1) 
+        Y = (k@X)/K[:,None]
+        return Y
+
+    #潜在変数の推定
+    def estimate_z(self, Z1, Z2, X, Y):
+        #ガウス関数による重み行列k(z_n,z_i)(NxN)
+        k = np.exp((-1/(2*(self.sigma**2)))*self.distance(Z1, Z2))
+        #kをiでsumした行列K(z_n)(Nx1)
+        K = np.sum(k, axis=1)
+
+        r = k/K[:,None]
+        d = Y[:,None,:] - X[None,:,:]
+        delta = Z1[:,None,:] - Z2[None,:,:]
+
+        B = Y - X
+        A = np.sum(B[:,None,:]*d,axis=2)
+
+        C = r*A
+        CC = C+C.T
+        E_bibun = (2/(self.N*(self.sigma**2)))*np.sum(CC[:,:,None]*delta,axis=1)
+
+        Z_new = Z1 - self.eta*(E_bibun + 2 * self.rambda * Z1 / self.N)
+
+        return Z_new
+
+    def fit(self, X, num_epoch, seed, resolution):
+        np.random.seed(seed)
+        Z = np.random.normal(scale=self.scale, size=(self.N, self.L))
+
+        history = dict(
+            E=np.zeros((num_epoch, )),
+            Y=np.zeros((num_epoch, self.N, self.D)),
+            f=np.zeros((num_epoch, resolution**self.L, self.D)),
+            Z=np.zeros((num_epoch, self.N, self.L))
+        )
+
+        with tqdm(range(num_epoch)) as pbar:
+                for epoch,i in enumerate(pbar):
+                    Y = self.estimate_y(Z, Z, X)
+                    Z = self.estimate_z(Z, Z, X, Y)
+
+                    """
+                    A = np.linspace(Z[:,0].min(), Z[:,0].max(), resolution)
+                    B = np.linspace(Z[:,1].min(), Z[:,1].max(), resolution)
+                    XX, YY = np.meshgrid(A,B)
+                    M = np.concatenate([XX.reshape(-1)[:,None], YY.reshape(-1)[:,None]], axis=1)
+                    """
+                    #sampleプログラムに合わせた初期値，make_gridの中身はよくわかってない
+                    Z_new = make_grid(resolution,
+                                    bounds=(np.min(Z), np.max(Z)),
+                                    dim=self.L)
+                    f = self.estimate_y(Z_new, Z, X)
+
+                    #E = (1/self.N) * np.sum(np.sum((Y - X)**2, axis=1) + self.rambda * np.sum(Z**2, axis=1), axis=0)
+                    E = (1/self.N) * (np.sum((Y - X)**2) + self.rambda * np.sum(Z**2))
+
+                    pbar.set_postfix(OrderedDict(loss=f'{E:.3f}'))
+
+                    history['E'][epoch] = E
+                    history['Y'][epoch] = Y
+                    history['f'][epoch] = f
+                    history['Z'][epoch] = Z
+
+        return history
+
+
+if __name__ == '__main__':
+    import data
+    from visualizer import visualize_history
+    N = 400
+    D = 3
+    L = 2
+    T = 100
+    eta = 10
+    sigma = 0.1
+    rambda = 0
+    scale=1e-8
+    resolution = 20
+
+    X = data.gen_saddle_shape(num_samples=N, random_seed=seed, noise_scale=0.05)
+    ukr = UKR(N, D, L, eta, sigma, rambda, scale)
+    history = ukr.fit(X, num_epoch=T, seed=seed, resolution=resolution)
+    visualize_history(X, history['Z'], history['f'], history['E'], color=X[:, 0], resolution=resolution, save_gif=True, 
+                      filename=f"saddle_N_{N}_D_{D}_L_{L}_T_{T}_eta_{eta}_sigma_{sigma}_scale_{scale}_resolution_{resolution}_seed_{seed}")