Multiple (not independent) response variables in machine learning Question: How to predict the percentage of people with age < 18, 18-65 and > 65 who visit a webpage using machine learning in R? Since these percentages sum to 100 for all observations, they are not independent.
Data: Contains 400 observations of 16 independent variables (TF-IDF values of 16 keywords for all pages) along with % of people with age < 18, 18-65 and > 65. This can be split into test and train set.
My idea: Recode the output variable as a factor (1: percentage of < 18 age group is highest, 2: percentage of 18-65 age group is highest, 3: otherwise). Build a generative classifier (like Naive Bayes) which predicts the probability of belonging to a each of the classes, which can be converted into percentage.
The problem: Variance of the output variable decreases on recoding. This may give inaccurate results.
 A: When you sum the probabilities of falling into each level of the outcome variable, of course they'll sum to 1, because the union of all the levels are your sample space, and P[sample space]=1. It has nothing to do with dependence, you actually can't even talk about dependence since your levels are not variables on their own, but just categories of the same target variable. In your case, the target is the percentage of visitors falling into each age group. In hydrology, you would have stream flow magnitude (low, medium, high), in marketing you would have sentimentality of a tweet (angry, sad, happy), and so forth.
Turning the output variable into a factor as you suggested makes sense, but it is true that it makes you lose some information (also, I think in that case you should NOT use class probabilities as proxies for percentages). To address this issue, why don't you keep your outcome variable continuous (pseudo-continuous actually since it can only take on discrete values I guess), and use a one-versus all approach where each ML algo learns to predict the percentage of a given age group based on all input variables. It would be a regression task rather than a classification task. Then, for a new observation, you can directly predict the percentage of each group, and aggregate the results.
One potential issue with this approach is that the percentages might exceed 100 in some instances. You could try to come up with some aggregation rule, or try to train a model to predict directly a 3 dimensional outcome (essentially you're trying to predict where to fall in the portion of a $100*100*100$ cube below the $x+y+z=100$ plane). It might be possible. I have found this reference about regression with a two-dimensional output variable which seems useful. The excellent lectures about spatial stats by Douglas Nychka from the NCAR might also be helpful. Some R code is also available.
A: Consider a binomial model instead of a model with several independent Bernoulli trials. Each row contains a summary: number of trials and number of successes for each combination of independent variables. Maximum likelihood logistic regression can be modified to take the sample size and number of successes as inputs.
The scenario given in the question is similar to the summarized form of binomial logistic regression where each input combination is repeated several times. Therefore, each row will have sample size and number of observations for each class of the multinomial distribution. Now multinomial logistic regression can be used as @EdM suggested. The predicted output will be a probability distribution over the set of possible outcome classes (age < 18, 18-65 and > 65). It is similar to applying multinomial logistic regression on a data set with observations repeated several times (total = sample size) with cumulative relative frequencies of outcomes matching the observed probabilities.
