View Single Post
  #173 (permalink)  
Old 05-13-2008, 07:59 AM
barnix barnix is offline
Senior Member
 
Join Date: Mar 2006
Posts: 716
barnix is on a distinguished road
using a Multi Layer Perceptron
Sébastien Marcel - Lab: Solving the Canon-ball problem using a Multi Layer Perceptron
Quote:
/* This program implements a simple Multi Layer Perceptron to solve the canonball problem (regression)

Adapted by Sebastien Marcel (2003-2004) from D. Collobert

The goal is to estimate the position X of a bullet at time t (fixed) given the initial speed (v) and angle (a)

y

/ \
|
|
|
| . X
| .
| .
| .
| .
| .
| .
| .
| .
| .
|. (v,a)
--------------------------> x

*/

#include <stdio.h>
#include <stdlib.h>
#include <math.h>

//
// Datasets
#define N_PATTERNS_TRAIN 500
#define N_PATTERNS_TEST 500

//*****************************
// Inputs
//
struct DataIn
{
float v; // speed of bullet
float a; // angle of bullet
} X[N_PATTERNS_TRAIN + N_PATTERNS_TEST];

//*****************************
// Targets
//
struct DataOut
{
float x; // x position at time t
float y; // y position at time t
} Y[N_PATTERNS_TRAIN + N_PATTERNS_TEST];

//*****************************
// to create the data
//
const float vmax = 1.0/10;
const float gravity = 1.0/200;

//*****************************
// to train the MLP
//
const float lambda = 0.01; // learning rate
const float mu = 0.6; // inertia momentum rate
const float mse_min = 0.001; // minimum Mean Squared Error
const int max_iterations = 500; // maximum number of iterations

//*****************************
// MLP data
//
#define N_HU 3

// weights between hidden neurons and inputs
float w1[N_HU + 1][3];
float w1old[N_HU + 1][3];

// weights between outputs and hidden neurons
float w2[2][N_HU + 1];
float w2old[2][N_HU + 1];

// values of integration function
float aHidden[N_HU + 1];
float aOutput[2];

// values of transfert function
float yHidden[N_HU + 1];

// constant value of bias
float xBias = 1.0;


//*****************************
// Sigmoid transfer function
//
float f(float x)
{
return (1.0 /(1.0 + exp(-x)));
}

//*****************************
// Derivative of Sigmoid
//
float f_prime(float x)
{
float z = exp(-x);
float one_plus_z = 1.0 + z;

return (z / (one_plus_z * one_plus_z));
}

#define MY_PI 3.141592

//*****************************
// Compute the x position of bullet at time t
//
float xpos(float t, float v, float teta)
{
return v * t * cos(teta*MY_PI/180.0);
}

//*****************************
// Compute the y position of bullet at time t
//
float ypos(float t, float v, float teta)
{
return -0.5 * gravity * t * t + v * t * sin(teta*MY_PI/180.0);
}

//*****************************
// Compute the MSE
//
float MSE(struct DataOut a, struct DataOut b)
{
return 0.5*((a.x - b.x)*(a.x - b.x) + (a.y - b.y)*(a.y - b.y));
}

//*****************************
// Random
//
float Random(float inf, float sup)
{
return( ((sup - inf)*((float)(random() & 0x7fffffffL) / (float) 0x7fffffffL) ) + inf);
}

//*****************************
// Create the MLP
//
void createMLP()
{
int i, j;
for (i = 0; i<=2; i++)
{
for (j = 0; j<= N_HU; j++)
{
w1[j][i] = Random(-1.0, 1.0);
w1old[j][i] = w1[j][i];
}
}

for (i = 0; i<2; i++)
{
for (j = 0; j<= N_HU; j++)
{
w2[i][j] = Random(-1.0, 1.0);
w2old[i][j] = w2[i][j];
}
}
}

//*****************************
// Create the datasets
//
void createDatasets()
{
int i;

i = 0;
while (i < (N_PATTERNS_TRAIN + N_PATTERNS_TEST))
{
X[i].v = (Random(0.0, vmax)) / vmax;
X[i].a = (Random(0.0, 90.0)) / (float)90;

Y[i].x = xpos(10, X[i].v * vmax, (X[i].a)*90.0);
Y[i].y = ypos(10, X[i].v * vmax, (X[i].a)*90.0);

if (Y[i].y > 0.0) i = i + 1;
}
}

//*****************************
// Forward the input in the MLP
//
void forward(struct DataIn In, struct DataOut *out)
{
int i;

aHidden[0] = xBias;
yHidden[0] = xBias;

for (i = 1; i <= N_HU; i++)
{
aHidden[i] = xBias*w1[i][0] + (In.v)*w1[i][1] + (In.a)*w1[i][2];
yHidden[i] = f(aHidden[i]);
}

aOutput[0] = 0.0;
aOutput[1] = 0.0;

for (i = 0; i <= N_HU; i++)
{
aOutput[0] = aOutput[0] + yHidden[i]*w2[0][i];
aOutput[1] = aOutput[1] + yHidden[i]*w2[1][i];
}

out->x = f(aOutput[0]);
out->y = f(aOutput[1]);
}

//*****************************
// Backward (back-propagate the gradient of error)
//
void backward(struct DataOut Yestimate, struct DataOut Outwant, struct DataIn In, float lambda, float mu)
{
float wnew, Goutput[2], Ghidden[N_HU + 1];
int i,j;

//
Goutput[0] = (Yestimate.x - Outwant.x) * f_prime(aOutput[0]);
Goutput[1] = (Yestimate.y - Outwant.y) * f_prime(aOutput[1]);

//
for (i = 0; i<=N_HU; i++)
Ghidden[i] = (Goutput[0]*w2[0][i] + Goutput[1]*w2[1][i])*f_prime(aHidden[i]);

//
for (j = 0; j<2; j++)
{
for (i = 0; i<=N_HU; i++)
{
wnew = w2[j][i] - lambda * yHidden[i] * Goutput[j] + mu * (w2[j][i] - w2old[j][i]);
w2old[j][i] = w2[j][i];
w2[j][i] = wnew;
}
}

//
for (i = 0; i<=N_HU; i++)
{
wnew = w1[i][0] - lambda * xBias * Ghidden[i] + mu * (w1[i][0] - w1old[i][0]);
w1old[i][0] = w1[i][0];
w1[i][0] = wnew;

wnew = w1[i][1] - lambda * In.v * Ghidden[i] + mu * (w1[i][1] - w1old[i][1]);
w1old[i][1] = w1[i][1];
w1[i][1] = wnew;

wnew = w1[i][2] - lambda * In.a * Ghidden[i] + mu * (w1[i][2] - w1old[i][2]);
w1old[i][2] = w1[i][2];
w1[i][2] = wnew;
}
}

//*****************************
// Stochastic gradient training
//
float train(int P)
{
struct DataOut Yestimate;
float mse_;
float mse_total;

mse_total = 0.0;

// For each train patterns
for(int p = 0 ; p < P ; p++)
{
// Forward current train pattern into the MLP
forward(X[p], &Yestimate);

// Computes the MSE
mse_ = MSE(Y[p], Yestimate);

// Accumulate the MSE
mse_total = mse_total + mse_;

// Backward MSE gradient
backward(Yestimate, Y[p], X[p], lambda, mu);
}

// Return normalized train MSE
return mse_total / (float) P;
}

//
int main()
{
float mse_train;
float mse_test;
float mse_;
struct DataOut Yestimate;

//*****************************
//
createMLP();

//*****************************
//
createDatasets();

//for (int i = 0; i< N_PATTERNS_TRAIN; i++) printf(" TRN: x=[%f %f] y=[%f %f]\n", X[i].v, X[i].a, Y[i].x, Y[i].y);
//for (int i = N_PATTERNS_TRAIN; i<(N_PATTERNS_TRAIN + N_PATTERNS_TEST); i++) printf(" TST: x=[%f %f] y=[%f %f]\n", X[i].v, X[i].a, Y[i].x, Y[i].y);

printf("Stochastic gradient training:\n");

//
FILE *pf_train = fopen("mse_train.txt", "w");
FILE *pf_test = fopen("mse_test.txt", "w");

int iter;

for(iter = 1 ; iter <= max_iterations ; iter++)
{
//*****************************
//
// train
mse_train = train(N_PATTERNS_TRAIN);

fprintf(pf_train, "%f\n", mse_train);


//*****************************
//
// test
mse_test = 0.0;

// For each test pattern
for(int i = N_PATTERNS_TRAIN ; i < (N_PATTERNS_TRAIN + N_PATTERNS_TEST) ; i++)
{
Yestimate.x = 0.0;
Yestimate.y = 0.0;

// forward current pattern into the MLP
forward(X[i], &Yestimate);

// computes MSE
mse_ = MSE(Y[i], Yestimate);

// accumulate MSE
mse_test += mse_;
}

// Normalize the MSE
mse_test /= (float)N_PATTERNS_TEST;

fprintf(pf_test, "%f\n", mse_test);
printf("."); fflush(stdout);

if (mse_train < mse_min) break;
}
printf("\n");

//
fclose(pf_test);
fclose(pf_train);

//*****************************
// Print info about training
//
printf("Number of iterations = %d\n", iter);
printf("Final MSE train = %f\n", mse_train);
printf("Final MSE test = %f\n", mse_test);

printf("End of program.\n");

return 0;
}
Attached Images
File Type: jpg canonball.jpg (5.6 KB, 813 views)

Last edited by barnix; 05-13-2008 at 08:01 AM.
Reply With Quote