What's the point in neural networks for multivariate regression? Do you have any case in which fitting a multivariate regression (so having multiple output nodes) outperforms the fitting of a single output one at a time in terms of accuracy? 
I ask this because there is not, as I see it, a sharing of information between output nodes, as could be the case in a linear multivariate regression where 
residuals of different responses can be correlated: 
$Y_{i1} = X_i \beta_1 + \epsilon_{i1}$ 
$Y_{i2} = X_i \beta_2 + \epsilon_{i2}$
$E[\epsilon_{i1}]= E[\epsilon_{i2}] = 0$, 
$Cov(\epsilon_{i1},\epsilon_{j2})=0$ for $i \neq j$, 
$Cov(\epsilon_{i1},\epsilon_{j2})= \sigma_{ij}$ for $i = j$ .  
In this simple (bivariate) case if $|\sigma_{ij}| \approx 1$ then I may expect that a very high value of $Y_1$ (relatively to its expected value given by $X \beta_1$) will result in a corresponding high value of $Y_2$ (relatively to its expected value $X \beta_2$ ). This is something that I don't find in a neural network context, since nodes are no random variables, rather deterministic.
 A: Here's a cartoon representation of the models. Model A uses a single network to predict both outputs. Model B uses separate networks to predict both outputs. Because the input to both networks is the same, we can re-frame this as the equivalent model C. In this case, the input is passed to to hidden/output layers that are simply concatenated copies of those in model B. The hidden layer weights would have a block diagonal structure, such that weights between units in the left and right halves are zero.

Presumably the output units are linear (because this is a regression problem) and the hidden units are nonlinear (otherwise why bother with a neural net). We can think of a network as mapping the input nonlinearly into a feature space. The images of the inputs in feature space are given by the activations of the last hidden layer. The output layer then performs linear regression in feature space. Training the network amounts to jointly learning the regression weights and feature space mapping.
Training model A would mean learning a single feature space mapping. Another way to say this is that we want to find a single representation of the input that's good for predicting both outputs. With model B, we'd find a separate representation for each output. By the equivalence to model C, this can also be seen as finding a single, higher dimensional representation.
Models B/C are much bigger networks than model A and have many more parameters. Consequently, they should be much more flexible. This could be good or bad depending on the situation. A bigger network can learn more complicated functions, given enough data. Given insufficient data, it can be more prone to overfit. Consider what would happen if we scaled this to 100 outputs instead of 2.
