python - Pythonとnumpyを使用した勾配降下

Question

def gradient(X_norm,y,theta,alpha,m,n,num_it):
    temp=np.array(np.zeros_like(theta,float))
    for i in range(0,num_it):
        h=np.dot(X_norm,theta)
        #temp[j]=theta[j]-(alpha/m)*(  np.sum( (h-y)*X_norm[:,j][np.newaxis,:] )  )
        temp[0]=theta[0]-(alpha/m)*(np.sum(h-y))
        temp[1]=theta[1]-(alpha/m)*(np.sum((h-y)*X_norm[:,1]))
        theta=temp
    return theta



X_norm,mean,std=featureScale(X)
#length of X (number of rows)
m=len(X)
X_norm=np.array([np.ones(m),X_norm])
n,m=np.shape(X_norm)
num_it=1500
alpha=0.01
theta=np.zeros(n,float)[:,np.newaxis]
X_norm=X_norm.transpose()
theta=gradient(X_norm,y,theta,alpha,m,n,num_it)
print theta

上記のコードからの私のシータは ですが、正しい matlab にある100.2 100.2はずです。100.2 61.09

Answer 1

あなたのコードは少し複雑すぎて、もっと構造が必要だと思います。最終的に、この回帰は 4 つの操作に要約されます。

1. 仮説 h = X * シータを計算します
2. 損失 = h - y を計算し、コストの 2 乗 (損失 ^2)/2m を計算します。
3. 勾配を計算する = X' * 損失 / m
4. パラメータを更新します theta = theta - alpha * gradient

mあなたの場合、 と混同していると思いますn。ここでmは、特徴の数ではなく、トレーニング セット内の例の数を示します。

あなたのコードの私のバリエーションを見てみましょう：

import numpy as np
import random

# m denotes the number of examples here, not the number of features
def gradientDescent(x, y, theta, alpha, m, numIterations):
    xTrans = x.transpose()
    for i in range(0, numIterations):
        hypothesis = np.dot(x, theta)
        loss = hypothesis - y
        # avg cost per example (the 2 in 2*m doesn't really matter here.
        # But to be consistent with the gradient, I include it)
        cost = np.sum(loss ** 2) / (2 * m)
        print("Iteration %d | Cost: %f" % (i, cost))
        # avg gradient per example
        gradient = np.dot(xTrans, loss) / m
        # update
        theta = theta - alpha * gradient
    return theta


def genData(numPoints, bias, variance):
    x = np.zeros(shape=(numPoints, 2))
    y = np.zeros(shape=numPoints)
    # basically a straight line
    for i in range(0, numPoints):
        # bias feature
        x[i][0] = 1
        x[i][1] = i
        # our target variable
        y[i] = (i + bias) + random.uniform(0, 1) * variance
    return x, y

# gen 100 points with a bias of 25 and 10 variance as a bit of noise
x, y = genData(100, 25, 10)
m, n = np.shape(x)
numIterations= 100000
alpha = 0.0005
theta = np.ones(n)
theta = gradientDescent(x, y, theta, alpha, m, numIterations)
print(theta)

最初に、次のような小さなランダム データセットを作成します。

ご覧のとおり、生成された回帰直線と、Excel で計算された数式も追加しました。

勾配降下法を使用した回帰の直感に注意する必要があります。データ X に対して完全なバッチ パスを実行する場合、すべての例の m 損失を 1 回の重み更新に減らす必要があります。この場合、これは勾配全体の合計の平均であるため、 で除算されmます。

次に注意する必要があるのは、収束を追跡し、学習率を調整することです。さらに言えば、反復ごとに常にコストを追跡し、プロットすることもできます。

私の例を実行すると、返される theta は次のようになります。

Iteration 99997 | Cost: 47883.706462
Iteration 99998 | Cost: 47883.706462
Iteration 99999 | Cost: 47883.706462
[ 29.25567368   1.01108458]

これは実際には、Excel で計算された式 (y = x + 30) に非常に近いものです。バイアスを最初の列に渡したので、最初の theta 値はバイアスの重みを示していることに注意してください。

Answer 2

以下に、線形回帰問題に対する勾配降下の実装を示します。

最初に、同様に勾配を計算X.T * (X * w - y) / Nし、この勾配で現在のシータを同時に更新します。

• X: 機能マトリックス
• y: 目標値
• w: 重み/値
• N: トレーニング セットのサイズ

Pythonコードは次のとおりです。

import pandas as pd
import numpy as np
from matplotlib import pyplot as plt
import random

def generateSample(N, variance=100):
    X = np.matrix(range(N)).T + 1
    Y = np.matrix([random.random() * variance + i * 10 + 900 for i in range(len(X))]).T
    return X, Y

def fitModel_gradient(x, y):
    N = len(x)
    w = np.zeros((x.shape[1], 1))
    eta = 0.0001

    maxIteration = 100000
    for i in range(maxIteration):
        error = x * w - y
        gradient = x.T * error / N
        w = w - eta * gradient
    return w

def plotModel(x, y, w):
    plt.plot(x[:,1], y, "x")
    plt.plot(x[:,1], x * w, "r-")
    plt.show()

def test(N, variance, modelFunction):
    X, Y = generateSample(N, variance)
    X = np.hstack([np.matrix(np.ones(len(X))).T, X])
    w = modelFunction(X, Y)
    plotModel(X, Y, w)


test(50, 600, fitModel_gradient)
test(50, 1000, fitModel_gradient)
test(100, 200, fitModel_gradient)