Where does the sum of squared errors function in neural networks come from? Training a basic multilayer perceptron neural network boils down to minimizing some kind of error function. Often the sum of squared errors is chosen as a this error function, but where does this function come from?

I always thought this function was chosen because it makes sense intuitively. However, recently I learned that this is only partly true and there is more behind it.

Bishop wrote in one of his papers that the sum of squared errors function can be derived from the principle of maximum likelihood. Furthermore he wrote that the squared error therefore makes the assumption that the noise on the target value has a Gaussian distribution.

I am not sure what he means with that. How does the sum of squared errors relate to the maximum likelihood principle in the context of neural networks?
you can trace the squared error in statistics through multivariate calculus all the way to Pythagorus. You are basically calculating the 'effective length' of the error, the hypotenuse, among errors from multiple variables $(X_1 - X_2)^2 + (Y_1 - Y_2)^2 + ...$ like in a triangle.

But where did the square root go?

Somebody realized that calculating roots of multiple variables over multiple iterations is computationally very expensive. So they decided to drop it. Checkout the squared Euclidean distance here for more details

How would a cubic error or a logarithmic error affect the outcome?

It just takes more time to converge because they are not as accurate. But we do see logarithmic errors over squares such as logistic regression where it is more optimal

All in all it is a simple case of optimization