@larsmans , Thanks for your observation on loss functions. 

Yes, the logistic loss  `-np.sum(Y * np.log(a_output) +  (1 - Y) * np.log(1 - a_output)) / a_output.shape[0]`
like you said, works for BOTH cases with expected accuracy, binary and multi class.  

The log loss `-log P(yt|yp) = -(yt log(yp) + (1 - yt) log(1 - yp))`, however, mine and the one implemented in sklearn, yield dire, different results on some binary-class problems. For logistic-based lbfgs on a binary-class toy example, my implementation achieves 0.8889, and sklearn's `log_loss` achieves 0.5. Where the actual score is in fact `1.000` when using `binary logistic loss`.

There must be an explanation for this discrepancy. To reproduce the scores, please see the file below, it contains two toy examples - binary and multi - with a test case. I run it whenever I change `mlp.py` lest I make wrong alterations.

https://dl.dropboxusercontent.com/u/30688032/mlp_example.py