NN.r 4（とりあえず完成版）

# calculate FeedForward Network
# Online version

# usage
# > source("NNonline.r")
# > init(nin, nhidden, nout, fact, fout) -> some global variables are defined here
# > bprop(x, y) -> learn by backpropagation
# > fprop(x) -> calculate feedforward propagation
#
# for debug
# > c <- test(hidden=10, times=2000, eta=0.01)
# > (plot of declining errors will be drawn)
# > n <- 1:length(c$y[,1])                     # means times
# > persp(n, c$x, c$y, theta=30, phi=30)

# names below are reserved
# (---variables---)
# dim.in, dim.hidden, dim.out, D, M, K, w1, w2, f.act, f.act.d, f.out
# (---functions---)
# init, fprop, bprop, sigmoid, identity

# functions
sigmoid <- function(x) 1/(1+exp(-x))
sigmoid.d <- function(x) sigmoid(x)*( 1 - sigmoid(x) )
identity <- function(x) x
identity.d <- function(x) return(1)

D <- numeric(0)
M <- numeric(0)
K <- numeric(0)
dim.in <- numeric(0)
dim.hidden <- numeric(0)
dim.out <- numeric(0)
w1 <- matrix(0)
w2 <- matrix(0)
f.act <- sigmoid
f.act.d <- sigmoid.d
f.out <- identity
f.out.d <- identity.d

# initializer
init <- function(nin, nhidden, nout, fact=sigmoid, fact.d=sigmoid.d, fout=identity, fout.d=identity.d){
# dimensions of input, hidden and output layers
D <<- dim.in <<- nin
M <<- dim.hidden <<- nhidden
K <<- dim.out <<- nout

f.act <<- fact
f.act.d <<- fact.d
f.out <<- fout
f.out.d <- fout.d

# initialize weights conforming to the PRML
w1 <<- matrix( runif( M*(D+1) ), M, D+1 ) # w1[,D+1] is bias
w2 <<- matrix( runif( K*(M+1) ), K, M+1 ) # w2[,M+1] is bias
}

# calculate forward propergation
# x : D dim. input matrix
fprop <- function( x ){
# calculate hidden layers
x.ex <- c(x, 1)
a1 <- w1 %*% x.ex
z <- f.act(a1)

# calculate output layers
z.ex <- c(z, 1)
a2 <- w2 %*% z.ex
y <- f.out(a2)

return( list( out=as.vector(y), out.a=as.vector(a2), hidden=as.vector(z), hidden.a=as.vector(a1)) )
}

# back propergation
# for sigmoid, identity only!
# calculate gradients of w1 and d2 for a data n
# data.in : input of training data
# data.out : output of training data
bprop <- function(data.in, data.out, eta=0.1){
# FP
calc <- fprop(data.in)

# calculate output layers -> renew w2
d2 <- calc$out - data.out
w2 <<- w2 - eta * ( d2 %o% c(calc$hidden, 1) )

# calculate hidden layers -> renew w1
d1 <- d2 %*% (w2[, -(M+1), drop=F]) %*% diag( f.act.d(calc$hidden.a) )
w1 <<- w1 - eta * ( as.vector(d1) %o% c(data.in, 1) )

if(any(is.na(w2))){
   print("there is NA in w2")
   return(list(d2, calc))
} else if(any(is.na(w1))) {
   print("there is NA in w1")
   return(list(d1, data.in))
}

return ( sqrt(sum(d2^2)) )
}

###########
# debug program

test <- function(times, hidden, eta){

init(1,hidden,1)
x <- seq(0, 3.14, by=0.1)
t <- sin(x)
e <- NULL
y <- NULL

for(i in 1:times){
  for(n in 1:length(x)){
    e.out <- bprop(x[n],t[n],eta)
    if(!is.numeric(e.out)){
      print(i,n)
      break
    }
     e <- c(e, e.out)
  }
  if(!is.numeric(e.out)){
    break
  }
  tmp <- NULL
  for(p in x){
   q <- fprop(p)
   tmp <- c(tmp, q$out)
  }
  y <- rbind(y, tmp)
}

plot(e)
return (list(error=e, x=x, t=t, y=y))
}