In [1]:

    

# まずはtheano tutorialの多項ロジット
class LogisticRegression(object):

    def __init__(self, input, n_in, n_out):

        self.W = theano.shared(
            value=numpy.zeros(
                (n_in, n_out),
                dtype=theano.config.floatX
            ),
            name='W',
            borrow=True
        )

        self.b = theano.shared(
            value=numpy.zeros(
                (n_out,),
                dtype=theano.config.floatX
            ),
            name='b',
            borrow=True
        )

        self.p_y_given_x = T.nnet.softmax(T.dot(input, self.W) + self.b)


        self.y_pred = T.argmax(self.p_y_given_x, axis=1)

        self.params = [self.W, self.b]

        self.input = input

    def negative_log_likelihood(self, y):

        # T.log(self.p_y_given_x)[T.arange(y.shape[0]), y]について
        # このスライスの仕方は下に例をのせる
        # ここでは、各サンプルに対して実際に得られているクラスへの判別確率を入手している
        return -T.mean(T.log(self.p_y_given_x)[T.arange(y.shape[0]), y])
       

    def errors(self, y):

        if y.ndim != self.y_pred.ndim:
            raise TypeError(
                'y should have the same shape as self.y_pred',
                ('y', y.type, 'y_pred', self.y_pred.type)
            )

        if y.dtype.startswith('int'):

            return T.mean(T.neq(self.y_pred, y))
        else:
            raise NotImplementedError()


    

In [18]:

    

# 二項ロジットで用いた関数をシグモイド関数と呼ぶ
# 多項ロジットで用いる関数をソフトマックス関数と呼ぶ
# 以下では人口データで多項ロジットを実装

% matplotlib inline
import numpy as np
import matplotlib.pyplot as plt

D = 2
N = 1500
K = 3

# データ点を用意
mean1 = [-2, 2]  # クラス1の平均
mean2 = [0, 0]   # クラス2の平均
mean3 = [2, -2]   # クラス3の平均
cov = [[1.0,0.8], [0.8,1.0]]  # 共分散行列（全クラス共通）

x1 = np.random.multivariate_normal(mean1, cov, N / 3)
x2 = np.random.multivariate_normal(mean2, cov, N / 3)
x3 = np.random.multivariate_normal(mean3, cov, N / 3)
x = np.vstack((x1, x2, x3))

# 教師ベクトルを用意
label1 = np.zeros((500,1))
label2 = np.zeros((500,1)) + 1
label3 = np.zeros((500,1)) + 2
labels = np.vstack((label1, label2, label3))


# 図示
plt.xlim((-6, 6))
plt.ylim((-6, 6))
plt.scatter(x1[:, 0], x1[:, 1], c='r')
plt.scatter(x2[:, 0], x2[:, 1], c='g')
plt.scatter(x3[:, 0], x3[:, 1], c='b')
plt.show()


    

In [19]:

    

# datasetにする
dataset = np.column_stack((x, labels))
np.random.shuffle(dataset) #データ点の順番をシャッフル


    

In [30]:

    

# theanoによる多項ロジット

import theano as th
import theano.tensor as T

# 定数項の追加
one = np.ones((1500,1))
x = np.column_stack((one, (dataset[:, :2])))
y = dataset[:, 2]
y = y.astype(np.int64)

tx = th.shared(np.asarray(x, dtype = th.config.floatX), borrow = True)
ty = th.shared(np.asarray(y, dtype = th.config.floatX), borrow = True)
theta = th.shared(np.zeros((3,3)), borrow = True)

# cost function
# hはサンプルごとに、各ラベルへの判別確率が出力されている行列
h = T.nnet.softmax(T.dot(tx, theta))
cost = -T.mean(T.log(h)[T.arange(y.shape[0]), y])

g_theta = T.grad(cost, theta)

learning_rate = 0.0001
updates = [(theta, theta - learning_rate * g_theta)]

train= th.function([], cost, updates = updates)

iteration = 100000
for itr in range(iteration):
    train()

print 'theta:', theta.get_value()


    

    




theta: [[-0.45296281  0.89566586 -0.44270305]
 [-0.97406508  0.00288094  0.97118414]
 [ 0.96121336  0.00909068 -0.97030404]]

In [13]:

    

b = np.array([[1,2,3,4,5],[6,7,8,9,10]])
a = np.array([2,3])
b[np.arange(2), a]


    

    
Out[13]:





array([3, 9])

In [22]:

    

T.arange(1500)


    

    
Out[22]:





ARange.0

In [26]:

    

y


    

    
Out[26]:





array([ 0.,  1.,  0., ...,  1.,  1.,  1.])

In [ ]: