Widorow-Hoffの学習規則

フリーソフトでつくる音声認識システム パターン認識・機械学習の初歩から対話システムまで

こちらの本のpp53-60にあるWidorow-Hoffの学習規則をRで実装しました。
いつも通り面倒なので理屈はほとんど省略
漢ならコードで語れ（笑）


基本的には前やったのとほぼ同じ。
パーセプトロンの学習規則は誤識別を許さないけど、今回のは誤識別を許す。
だから、線形分離不可能なデータでもある程度分類できる。
方法としては以下の式、つまり二乗和誤差関数の最小化を行う。

ただし、wをクラスの数だけ用意しないといけない。
2クラス分類(C=2)のときは２つ用意しないといけない。
wは下の更新式で求めていく。（最急降下法）

今回は前回と違ってベクトルでの表現が簡単だからRとの相性も良さそう。
はとすれば簡単にベクトル表現ができる。
以下Rでの実装コード

# Widrow-Hoffの学習規則
# 2乗和誤差の最小化
# 最小化には最急降下法を用いる
# 2クラス判別用
#識別に用いる関数
# x:入力ベクトル, w:重みベクトル
myWH.f <- function(x,w){
  x <- c(1,x)
  return ( sum(w*x) )
}
#クラス識別関数
# x:入力ベクトル, w1:重みベクトル1, w2:重みベクトル2
myWH.predict <- function(x,w1,w2){
  cls1 <- myWH.f(x,w1)
  cls2 <- myWH.f(x,w2)
  if( cls1 > cls2 ){
    return (1)
  }else{
    return (2)
  }
}
#学習関数
# gc(x)の係数wを求める(c:クラス)
# data:学習データ, t:学習データのクラスベクトル(1or2) )
# its:最大学習回数, c:識別するクラス
myWH.train <- function(data, t, its, c){
  # 学習の準備
  RHO.MIN <- 1e-8
  if( is.vector(data) == TRUE ){
    data <- matrix(data,length(data),1) # ベクトルを行列に直す
  }
  # 一方のクラスのクラスベクトルを0とする(計算のため)
  if( c == 1 ){#クラス2 のクラスベクトルを0とする
    t <- replace(t, t==2, 0)
  }else{       #クラス1のクラスベクトルを0とし、クラス2のクラスベクトルを1とする
    t <- replace(t, t==1, 0)
    t <- replace(t, t==2, 1)
  }
  # 学習データに0次元目(=1)を加える(biasの計算)
  data2 <- matrix( 1, length(data[,1]), length(data[1,])+1)
  for( i in 1:length(data[,1]) ){
    for( j in 1:length(data[1,]) ){
      data2[i,j+1] = data[i,j]
    }
  }
  data <- data2 
  #2乗和誤差関数
  # w:重みベクトル
  myWH.gosa <- function(w){
    yosoku <- function(x){
      return ( sum(w*x) )
    }
    gosa <- apply(data,1,yosoku) - t
    J = sum(gosa^2) # 1/2は省略
    return (J)
  }
  #勾配ベクトルを求める関数
  # w:重みベクトル(バイアスを含む)
  myWH.koubai <- function(w){
    yosoku <- function(x){
      return ( sum(w*x) )
    }
    gosa <- apply(data,1,yosoku) - t
    JJ = gosa%*%data
    if( is.matrix(JJ) ){
      JJ <- JJ[1,] #行列になっているのでベクトルに直す
    }
    return (JJ)
  }
  #ここから学習の本番
  w <- numeric(ncol(data))
  gosa.old <- myWH.gosa(w)
  cat(" 初期誤差=",gosa.old,"\n")
  for( i in 1:its ){
    rho <- 1
    while( rho > RHO.MIN ){
      JJ <- myWH.koubai(w)
      w.new <- (w) - (rho * JJ)
      gosa.new <- myWH.gosa(w.new)
      if( gosa.new > gosa.old ){
        rho <- rho / 2
      }else{
        w <- w.new
        break
      }
    }
    # 誤差の修正幅が閾値未満になると終了(収束条件)
    if( abs(gosa.old-gosa.new) < 1e-8 ){
      break
    }
    gosa.old <- gosa.new
  }
  gosa.final <- myWH.gosa(w)
  cat("収束回数 ",i,"回\n")
  cat("最終誤差 ",gosa.final,"\n")
  return (w)
}
#プロット関数 2次元2クラス用
# w1:学習結果1,w2:学習結果2, data:学習データ, t:学習データのクラスラベル
myWH.plot <- function(w1,w2,data,t){
  a <- 200 # プロットデータの個数,大きくすると時間がかかる
  tx <- seq(0,9,length=a)
  ty <- tx 
  rapMyWH <- function(x,y){
    v <- c(x,y)
    return ( myWH.predict(v,w1,w2) )
  }
  z <- matrix(0,length(ty),length(tx))
  for( y in 1:length(tx) ){
    for( x in 1:length(ty) ){
      z[x,y] <- rapMyWH(tx[x],ty[y])
    }
  }
  # 分類結果の領域を描画
  image(tx,ty,z,col=c("#00FF80","#FFBF00"),axes=TRUE,xlab="x",ylab="y")
  # 学習データ点をプロット
  points( datasrc[,1], datasrc[,2], col=ifelse(t==1,"red","blue"), pch=ifelse(t==1,16,17) )
  title("myWH:[testData0.csv]")
  return (z)
}

前回同様データはtestData0.csvを用いた。
実行コード

datasrc <- read.csv("testData0.csv")
t <- datasrc[,3]
datasrc <- datasrc[,1:2]
data.matrix <- matrix(0,nrow(datasrc),ncol(datasrc))
for( i in 1:nrow(datasrc) ){
  for( j in 1:ncol(datasrc) ){
    data.matrix[i,j] <- datasrc[i,j]
  }
}
its <- 1000 #最大学習回数
(w1 <- myWH.train(data.matrix, t, its, 1) )
(w2 <- myWH.train(data.matrix, t, its, 2) )
z <- myWH.plot(w1,w2,data,t) #学習結果の表示

実行した結果のプロット図

とりあえずきれいに分類できてるから良しとする。