How to determine the size of biases?
The biases appear in the activation function and can be seen as node feature. Suppose you have act = sigmoid, X your input layer and W your first hidden layer weight matrix.
Then the output matrix of this hidden layer will be act(X.dot(W)+b). 

So your bias should be a vector containing a random value for each node of the current hidden layer. 

It will shift randomly the weight and so the activation function. See here (usually it prevent from overfitting).

I didn't read the whole code so I don't know if the problem comes from here, but why don't you learn how to use tensorflow instead of re-write everything ? There are a lot of tutorials in the net, like MNIST.