For example let's say we want to simultaneously predict GDP and Population as a linear function of the Area and Coastline Length for each country having sufficient data:
data = Cases[ Table[CountryData[country, property], {country, CountryData[]}, {property, {"Area", "CoastlineLength", "GDP", "Population"}}], {_Real ..}];
There is no problem predicting either GDP or Population separately. Each of the following works:
LinearModelFit[data[[All, {1, 2, 3}]], {area, coast}, {area, coast}]
LinearModelFit[data[[All, {1, 2, 4}]], {area, coast}, {area, coast}]
And one would imagine this would generalize to multivariate regressions with multiple outcome variables naturally:
LinearModelFit[{#1, #2, {#3, #4}} & @@@ data, {area, coast}, {area, coast}]
But this comes back with an error:
LinearModelFit::notdata: The first argument is not a vector, matrix, or a list containing a design matrix and response vector. >>