NN.r 3

# 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

# 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)

print(list(w1,w2))

# 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) )
}

###########
# debug program
init(1,4,1)
x <- seq(0,3.14, by=0.1)
y <- sin(x)
for(i in 1:100){
for(n in 1:length(x))
  bprop(x[i],y[i])
}

p.y <- NULL
for(p in x){
q <- fprop(p)
p.y <- c(p.y, q$out)
}

plot(x,p.y)