if(p < 0 || p > 1) return 0;

int c = 0;

double r;

for(int i=0; i<n; i++) {

r = rand() / (RAND_MAX + 1.0);

if (r < p) c++;

}

return c;

N = size;

n_visible = n_v;

n_hidden = n_h;

if(w == NULL) {

W = new double*[n_hidden];

for(int i=0; i<n_hidden; i++) W[i] = new double[n_visible];

double a = 1.0 / n_visible;

for(int i=0; i<n_hidden; i++) {

for(int j=0; j<n_visible; j++) {

W[i][j] = uniform(-a, a);

}

}

} else {

W = w;

}

if(hb == NULL) {

hbias = new double[n_hidden];

for(int i=0; i<n_hidden; i++) hbias[i] = 0;

} else {

hbias = hb;

}

if(vb == NULL) {

vbias = new double[n_visible];

for(int i=0; i<n_visible; i++) vbias[i] = 0;

} else {

vbias = vb;

}

for(int i=0; i<n_visible; i++) {

if(x[i] == 0) {

tilde_x[i] = 0;

} else {

tilde_x[i] = binomial(1, p);

}

}

for(int i=0; i<n_hidden; i++) {

y[i] = 0;

for(int j=0; j<n_visible; j++) {

y[i] += W[i][j] * x[j];

}

y[i] += hbias[i];

y[i] = sigmoid(y[i]);

}

for(int i=0; i<n_visible; i++) {

z[i] = 0;

for(int j=0; j<n_hidden; j++) {

z[i] += W[j][i] * y[j];

}

z[i] += vbias[i];

z[i] = sigmoid(z[i]);

}

int *tilde_x = new int[n_visible];

double *y = new double[n_hidden];

double *z = new double[n_visible];

double *L_vbias = new double[n_visible];

double *L_hbias = new double[n_hidden];

double p = 1 - corruption_level;

get_corrupted_input(x, tilde_x, p);

get_hidden_values(tilde_x, y);

get_reconstructed_input(y, z);

// vbias

for(int i=0; i<n_visible; i++) {

L_vbias[i] = x[i] - z[i];

vbias[i] += lr * L_vbias[i] / N;

}

// hbias

for(int i=0; i<n_hidden; i++) {

L_hbias[i] = 0;

for(int j=0; j<n_visible; j++) {

L_hbias[i] += W[i][j] * L_vbias[j];

}

L_hbias[i] *= y[i] * (1 - y[i]);

hbias[i] += lr * L_hbias[i] / N;

}

// W

for(int i=0; i<n_hidden; i++) {

for(int j=0; j<n_visible; j++) {

W[i][j] += lr * (L_hbias[i] * tilde_x[j] + L_vbias[j] * y[i]) / N;

}

}

delete[] L_hbias;

delete[] L_vbias;

delete[] z;

delete[] y;

delete[] tilde_x;

double *y = new double[n_hidden];

get_hidden_values(x, y);

get_reconstructed_input(y, z);

delete[] y;

srand(0);

double learning_rate = 0.1;

double corruption_level = 0.3;

int training_epochs = 100;

int train_N = 10;

int test_N = 2;

int n_visible = 20;

int n_hidden = 5;

// training data

int train_X[10][20] = {

{1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},

{1, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},

{1, 1, 0, 1, 1, 1, 1, 1, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},

{1, 1, 1, 1, 0, 0, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},

{0, 1, 1, 1, 1, 1, 1, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},

{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1},

{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1},

{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 1, 1, 1, 1, 1, 0, 1},

{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 0, 0, 1, 1, 1, 1},

{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1, 1, 1, 0}

};

// construct dA

dA da(train_N, n_visible, n_hidden, NULL, NULL, NULL);

// train

for(int epoch=0; epoch<training_epochs; epoch++) {

for(int i=0; i<train_N; i++) {

da.train(train_X[i], learning_rate, corruption_level);

}

}

// test data

int test_X[2][20] = {

{1, 0, 1, 0, 1, 0, 1, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0},

{0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 1, 1, 0, 1, 1, 0}

};

double reconstructed_X[2][20];

// test

for(int i=0; i<test_N; i++) {

da.reconstruct(test_X[i], reconstructed_X[i]);

for(int j=0; j<n_visible; j++) {

printf("%.5f ", reconstructed_X[i][j]);

}

cout << endl;

}

cout << endl;