building a multilayer perceptron

[E093605o]

I want to build a multilayer perceptron and am almost done. I cannot initialize a test object though. I have written a basic linear algebra library "Matrix.h" to use.
In the main() of the MLP.cpp file I commented the error above the line producing it.

//MLP.h
#pragma once
#include "Matrix.h"

template<typename T>
class MLP {
 public:
  std::vector<size_t> units_per_layer;
  std::vector<Matrix<T>> bias_vectors;
  std::vector<Matrix<T>> weight_matrices;
  std::vector<Matrix<T>> activations;

  double lr = .001d;

  MLP(std::vector<size_t> units_per_layer):
    units_per_layer(units_per_layer),
    weight_matrices(),
    bias_vectors(),
    activations()
    {

  for (size_t i = 0; i < units_per_layer.size() - 1; ++i) {
    size_t in_channels{units_per_layer[i]};
    size_t out_channels{units_per_layer[i+1]};

    // initialize to random Gaussian
    auto W  = mtx<T>::randn(out_channels, in_channels);
    weight_matrices.push_back(W);

    auto b  = mtx<T>::randn(out_channels, 1);
    bias_vectors.push_back(b);

    activations.resize(units_per_layer.size());
  }
}

inline auto sigmoid(double x) {
  return 1.0 / (1 + exp(-x));
}

inline auto d_sigmoid(double x){
  return (x * (1 - x));
}   

template <typename T>
auto MLP::forward(Matrix<T> x) {
  assert(std::get<0>(x.shape) == units_per_layer[0] && std:::get<1>(x.shape));

  activations[0] = x;
  Matrix prev(x);
  for (int i = 0; i < units_per_layer.size() - 1; ++i) {

    Matrix y = weight_matrices[i].matmul(prev);
    y = y + bias_vectors[i];
    y = y.apply_function(sigmoid);
    activations[i+1] = y;
    prev = y;
  }
  return prev;
}

template<typename T>
void MLP<T>::backprop(Matrix<T> target) {
  assert(get<0>(target.shape) == units_per_layer.back());

  // determine the simple error
  // error = target - output
  auto y = target;
  auto y_hat = activations.back();
  auto error = (target - y_hat);

  // backprop the error from output to input and step the weights
  for(int i = weight_matrices.size() - 1 ; i >= 0; --i) {
    //calculating errors for previous layer
    auto Wt = weight_matrices[i].T();
    auto prev_errors = Wt.matmul(delta);

    // apply derivative of function evaluated at activations
    //backprop for biases
    auto d_outputs = activations[i+1].apply_function(d_sigmoid);
    auto gradients = error.multiply_elementwise(d_outputs);
    gradients = gradients.multiply_scalar(lr);

    // backprop for weights
    auto a_trans = activations[i].T();
    auto weight_gradients = gradients.matmul(a_trans);

    //adjust weights
    bias_vectors[i] = bias_vectors[i].add(gradients);
    weight_matrices[i] = weight_matrices[i].add(weight_gradients);
    error = prev_errors;
  }
}

};

//MLP.cpp

#include "Matrix.h"
#include "MLP.h" 


int main(){
    
    std::vector<size_t> layers = {3,3};
    /* the following line throws an error: no instance of constructor "MLP" matches 
    the argument list */
    MLP test(layers);
    return 0;
}

[kbw]

Can't see what the specific problem is, and I don't have your Matrix.h, I can't compile it.

However, symbol units_per_layer in the constructor is ambiguous; not for the compiler, but for you.

In the initializer list, you really should use:

  units_per_layer(std::move(units_per_layer))

as there's no reason to copy that vector. You should also move it in main to remove that copy:

  MLP test(std::move(layers));

When you refer to units_per_layer in the constructor, you're picking up the parameter, not the member. It's best to give the parameter a different name.

Distilling your code to:

#include <vector>

template<typename T>
class MLP {
        std::vector<size_t> units_per_layer;

public:
        MLP(std::vector<size_t> units_per_layer):
                units_per_layer(units_per_layer) {
        }
};

int main() {
        MLP<int>{{2, 3, 4}};
}

It's seems the problem is MLP is a template class, and you didn't specify the type.

[E093605o]

thanks for your answer. That indeed was one of the problems, it still does not compile however. The compiler complains "error: extra qualification ‘MLP<T>::’ on member ‘forward’ [-fpermissive]
45 | auto MLP<T>::forward(Matrix<T> x)" and the same for the other template method. Do you know why this happens? And the same error is raised for the backprop method as it has the same syntactical structure.

[seeplus]

There's some issues with templates in MLP. There's a different syntax for defining a member function separate to the class definition (Matrix) to defining within the class definition (MLP). Also, delta isn't defined. Also note that there are some copies which possibly aren't needed which could reduce performance - but I'd suggest getting the code to work first and only if the performance is not acceptable for your real world data to then look at further optimisation. I've changed the params for a couple of functions to be pass by const ref as opposed to by value to remove some copying of data data structures.

This compiles:

#include <vector>
#include <cmath>
#include <cassert>
#include <iostream>
#include <tuple>
#include <random>
#include <functional>

template<typename Type>
class Matrix {

	size_t cols {};
	size_t rows {};

public:
	std::vector<std::vector<Type>> data;
	std::tuple<size_t, size_t> shape;
	size_t elementCount {};

	/* constructors */
	Matrix(size_t rowsArg, size_t colsArg) : cols(colsArg), rows(rowsArg),
		elementCount(rows* cols), shape(std::tuple<size_t, size_t>(rows, cols)) {
		data = std::vector<std::vector<Type>>(rows, std::vector<Type>(cols));
	}

	Matrix() {};

	//methods
	void print();

	Matrix<Type> matmul(Matrix<Type>& m);

	Matrix<Type> multiply_elementwise(Matrix<Type>& m);

	Matrix<Type> multiply_scalar(Type scalar);

	Matrix<Type> square();

	Matrix<Type> add(Matrix<Type>& m);

	Matrix<Type> sub(Matrix& target);

	Matrix<Type> T();

	Matrix<Type> apply_function(Type(*func)(Type));

	Type& operator()(size_t row, size_t col) {
		assert(row < data.size() && col < data[0].size());
		return data[row][col];
	}

	Matrix operator+(Matrix& target) {
		return add(target);
	}

	Matrix operator-() {
		Matrix output(rows, cols);
		for (size_t r = 0; r < rows; ++r) {
			for (size_t c = 0; c < cols; ++c) {
				output(r, c) = -(*this)(r, c);
			}
		}
		return output;
	}

	Matrix operator-(Matrix& target) {  // for cleaner usage
		return sub(target);
	}


};

// methods
template<typename Type>
void Matrix<Type>::print() {
	for (int i = 0; i < rows; i++) {
		for (int j = 0; j < cols; j++) {
			std::cout << data[i][j] << " ";
		}
		std::cout << std::endl;
	}
}

template <typename Type>
Matrix<Type> Matrix<Type>::matmul(Matrix<Type>& target) {
	assert(cols == target.rows);
	Matrix output(rows, target.cols);

	for (size_t r = 0; r < output.rows; ++r) {
		for (size_t c = 0; c < output.cols; ++c) {
			for (size_t k = 0; k < target.rows; ++k)
				output(r, c) += (*this)(r, k) * target(k, c);
		}
	}
	return output;
};

template <typename T>
struct mtx {
	static Matrix<T> randn(size_t rows, size_t cols) {
		Matrix<T> M(rows, cols);

		std::random_device rd {};
		std::mt19937 gen { rd() };

		// init Gaussian distr. w/ N(mean=0, stdev=1/sqrt(numel))
		T n = static_cast<T>(M.elementCount);
		T stdev { 1 / sqrt(n) };
		std::normal_distribution<T> d { 0, stdev };

		// fill each element w/ draw from distribution
		for (size_t r = 0; r < rows; ++r) {
			for (int c = 0; c < cols; ++c) {
				M(r, c) = d(gen);
			}
		}
		return M;
	}
};

template <typename Type>
Matrix<Type> Matrix<Type>::multiply_elementwise(Matrix<Type>& target) {
	assert(shape == target.shape);
	Matrix output((*this));
	for (size_t r = 0; r < output.rows; ++r) {
		for (size_t c = 0; c < output.cols; ++c) {
			output(r, c) = target(r, c) * (*this)(r, c);
		}
	}
	return output;
}

template<typename Type>
Matrix<Type> Matrix<Type>::square() {
	Matrix output((*this));
	output = multiply_elementwise(output);
	return output;
}

template<typename Type>
Matrix<Type> Matrix<Type>::multiply_scalar(Type scalar) {
	Matrix output((*this));
	for (size_t r = 0; r < output.rows; ++r) {
		for (size_t c = 0; c < output.cols; ++c) {
			output(r, c) = scalar * (*this)(r, c);
		}
	}
	return output;
}

template<typename Type>
Matrix<Type> Matrix<Type>::add(Matrix& target) {
	assert(shape == target.shape);
	Matrix output(rows, std::get<1>(target.shape));

	for (size_t r = 0; r < output.rows; ++r) {
		for (size_t c = 0; c < output.cols; ++c) {
			output(r, c) = (*this)(r, c) + target(r, c);
		}
	}
	return output;
}

template<typename Type>
Matrix<Type> Matrix<Type>::sub(Matrix& target) {
	Matrix neg_target = -target;
	return add(neg_target);
}

template<typename Type>
Matrix<Type> Matrix<Type>::T() {
	size_t new_rows { cols }, new_cols { rows };
	Matrix transposed(new_rows, new_cols);

	for (size_t r = 0; r < new_rows; ++r) {
		for (size_t c = 0; c < new_cols; ++c) {
			transposed(r, c) = (*this)(c, r);  // swap row and col
		}
	}
	return transposed;
}

template<typename Type>
Matrix<Type> Matrix<Type>::apply_function(Type(*func)(Type)) {
	Matrix output((*this));
	for (size_t r = 0; r < rows; ++r) {
		for (size_t c = 0; c < cols; ++c) {
			output(r, c) = func((*this)(r, c));
		}
	}
	return output;
}


template<typename T>
class MLP {
public:
	std::vector<size_t> units_per_layer;
	std::vector<Matrix<T>> bias_vectors;
	std::vector<Matrix<T>> weight_matrices;
	std::vector<Matrix<T>> activations;

	double lr = .001;

	MLP(const std::vector<size_t>& units_per_layer) :
		units_per_layer(units_per_layer) {
		//weight_matrices(),
		//bias_vectors(),
		//activations() {

		for (size_t i = 0; i < units_per_layer.size() - 1; ++i) {
			size_t in_channels { units_per_layer[i] };
			size_t out_channels { units_per_layer[i + 1] };

			// initialize to random Gaussian
			auto W = mtx<T>::randn(out_channels, in_channels);
			weight_matrices.push_back(W);

			auto b = mtx<T>::randn(out_channels, 1);
			bias_vectors.push_back(b);

			activations.resize(units_per_layer.size());
		}
	}

	inline auto sigmoid(double x) {
		return 1.0 / (1 + exp(-x));
	}

	inline auto d_sigmoid(double x) {
		return (x * (1 - x));
	}

	auto forward(const Matrix<T>& x) {
		assert(std::get<0>(x.shape) == units_per_layer[0] && std:::get<1>(x.shape));

		activations[0] = x;
		auto prev(x);

		for (int i = 0; i < units_per_layer.size() - 1; ++i) {
			auto y = weight_matrices[i].matmul(prev);

			y = y + bias_vectors[i];
			y = y.apply_function(sigmoid);
			activations[i + 1] = y;
			prev = y;
		}
		return prev;
	}

	void backprop(const Matrix<T>& target) {
		assert(get<0>(target.shape) == units_per_layer.back());

		// determine the simple error
		// error = target - output
		auto y = target;
		auto y_hat = activations.back();
		auto error = (target - y_hat);

		// backprop the error from output to input and step the weights

		for (int i = weight_matrices.size() - 1; i >= 0; --i) {
		  //calculating errors for previous layer
			auto Wt = weight_matrices[i].T();

			//// delta NOT DEFINED
			//auto prev_errors = Wt.matmul(delta);

			// apply derivative of function evaluated at activations
			//backprop for biases
			auto d_outputs = activations[i + 1].apply_function(d_sigmoid);
			auto gradients = error.multiply_elementwise(d_outputs);
			gradients = gradients.multiply_scalar(lr);

			// backprop for weights
			auto a_trans = activations[i].T();
			auto weight_gradients = gradients.matmul(a_trans);

			//adjust weights
			bias_vectors[i] = bias_vectors[i].add(gradients);
			weight_matrices[i] = weight_matrices[i].add(weight_gradients);
			////error = prev_errors;
		}
	}
};


int main() {
	const std::vector<size_t> layers { 3,3 };

	MLP<double> test (layers);
}

[kbw]

Your indentation seems to be out, making it difficult to see that you've written:

template<typename Type>
class Matrix {
public:
    template <typename T>
    auto MLP::forward(Matrix<T> x) {
        //...
    }
};

The error should now be obvious. (That MPL:: shouldn't be there.)

I suggest you use tabs or a shiftwidth of 4, and ensure your indentation is correct. It can be checked with clang-tidy.

[E093605o]

thanks for all your answers. The MLP works now. One last thing, I want to log the error and outputs in a file to plot some graphs later, how do I do the logging part? I want to use fstream in the Main.cpp file but I did not find helpful information right away.

//MLP.h
#pragma once
#include "Matrix.h"

template<typename T>
class MLP {
 public:
  std::vector<size_t> units_per_layer;
  std::vector<Matrix<T>> bias_vectors;
  std::vector<Matrix<T>> weight_matrices;
  std::vector<Matrix<T>> activations;
  std::vector<Matrix<T>> zs;

  double lr = .001d;

  MLP(std::vector<size_t> units_per_layer):
    units_per_layer(units_per_layer),
    weight_matrices(),
    bias_vectors(),
    zs(),
    activations()
    {

  for (size_t i = 0; i < units_per_layer.size() - 1; ++i) {
    size_t in_channels{units_per_layer[i]};
    size_t out_channels{units_per_layer[i+1]};

    // initialize to random Gaussian
    auto W  = mtx<T>::randn(out_channels, in_channels);
    weight_matrices.push_back(W);

    auto b  = mtx<T>::randn(out_channels, 1);
    bias_vectors.push_back(b);

    auto z = mtx<T>::randn(out_channels,1);
    zs.push_back(z);

    activations.resize(units_per_layer.size());
  }
}

inline auto sigmoid(double x) {
  return 1.0 / (1 + exp(-x));
}

inline auto d_sigmoid(double x){
  return (x * (1 - x));
}   


auto forward(Matrix<T> x) {
  assert(std::get<0>(x.shape) == units_per_layer[0] && std::get<1>(x.shape));

  activations[0] = x;
  Matrix prev(x);
  for (int i = 0; i < units_per_layer.size() - 1; ++i) {

    Matrix y = weight_matrices[i].matmul(prev);
    y = y + bias_vectors[i];
    y = y.apply_function(sigmoid);
    activations[i+1] = y;
    prev = y;
  }
  return prev;
}

void backprop(Matrix<T> &target) {
  assert(std::get<0>(target.shape) == units_per_layer.back());

  // determine the simple error
  // error = target - output
  Matrix<T> y = target;
  Matrix<T> y_hat = activations.back();
  Matrix<T> error = (target - y_hat);
  Matrix<T> last_z = zs[zs.size()-1];
  Matrix<T> delta_L = error.multiply_elementwise(last_z.apply_function(d_sigmoid));

  // backprop the error from output to input and step the weights
  for(int i = weight_matrices.size() - 2 ; i >= 0; --i) {
    // calculating error for previous layer
    Matrix<T> delta_l = weight_matrices[i].T().matmul(delta_L).multiply_elementwise(zs[i].apply_function(d_sigmoid));

    // calculating the change of weights and biases
    Matrix<T> nabla_w = delta_l.matmul(activations[i].T());
    Matrix<T> nabla_b = delta_l;

    // updating weights and biases
    weight_matrices[i] = weight_matrices[i] - nabla_w.multiply_scalar(lr);
    bias_vectors[i] = bias_vectors[i] - delta_l.multiply_scalar(lr);

    // updating the error term delta
    delta_L = delta_l;
  }
}

};

//Main.cpp
#pragma once
#include "Matrix.h"
#include "MLP.h"
#include <vector>
#include <iostream>
#include <fstream>
#include <math.h>

int main() {

  // init model
  std::vector<size_t> layers = {1,8,8,8,1};


  // open file to save loss, x, y, and model(x)
  std::ofstream my_file; 
  my_file.open ("data.txt");

  int max_iter{10000};

  const double PI {3.14159};
  

  MLP<double> model(layers);

  for (int i = 0; i < max_iter; i++){
    auto x = mtx<double>::randn(1, 1).multiply_scalar(PI);
    auto y = x.apply_function([](double v) -> double { return sin(v) * sin(v); });

  // forward and backward
  auto y_hat = model.forward(x); 
  model.backprop(y); // loss and grads computed in here

  // function that logs (loss, x, y, y_hat)
  log(my_file, x, y, y_hat); 
  }

  my_file.close();
}

[kbw]

I want to log the error and outputs in a file to plot some graphs later, how do I do the logging part?

I'm not sure what you mean, just open an output file and write to it. You can write to std::clog for logging, that output will go to stderr. It's up to you, I'm not sure what you want to do.

However personally, I'd use spdlog for logging, and use a std::ofstream to write the data to a file.
https://github.com/gabime/spdlog

[E093605o]

the thing is, in line 133 it says "no instance of overloaded function "log" matches the argument list". Might it be the case that to log this file I need some function to have been implemented in Matrix.h?

[seeplus]

If I'm reading the code right x, y, y_hat are all of type Matrix?

So you need a function:

template<typename T>
void log (std::ostream& file, const Matrix<T>& x, const Matrix<T>& y, const Matrix<T>& y_hat) {
    // use file as the stream to write data to
}

Note L133 should be my_file and not file

eg like the Matrix print function but use file instead of std::cout

For each Matrix, how do you want the data written?

What program is going to use this data? Does it have a format requirement?

[E093605o]

Yes they're all of type Matrix. Once I logged this I thought of writing a small python script that plots these values to see if the neural net actually works.
I think it would be easiest to have one line per log.
Each Matrix denotes a vector here, input vector, true output vector and output vector of the NN. I would like them to be logged in one line.

[seeplus]

If you want one line per Matrix - and Matrix is a 'vector of vector', how do you want to separate each row? or start line with the number of row/cols or ??

When devising log formats, this often depends upon what's easiest for the reader. As I don't use python, what's the easiest format of logs for python to read?

[seeplus]

You might find it useful to overload the operator<< for the Matrix class eg:

template<typename T>
std::ostream& operator<<(std::ostream& os, const Matrix<T>& mat) {
    for (const auto& r : data) {
        for (const auto& e : r)
            os << e << ' ';

        os << '\n';
    }

    return os;
}

As data is public, this doesn't need to be a friend function.

Then log becomes something like:

template<typename T>
void log (std::ostream& file, const Matrix<T>& x, const Matrix<T>& y, const Matrix<T>& y_hat) {
    file << x << '\n';
    file << y << '\n';
    file << y_hat << '\n';
}

[E093605o]

sorry I went for a run, thats why I'm replying so late. I wanted to put all the contents of one log() call into one line, maybe thats easier to parse because you can then read contents line by line and strip special separation characters. I think the easiest format for python would be \n separated strings.