What are the appropriate Statistical tests when dealing with two or more categorical outcome variables in a model I'm familiar with GLM techniques when dealing with continuous outcomes and categorical or continuous predictors.
I've used Logistic regression with a single dichotomous categorical outcome variable and continuous predictor variables.
I thought today, what about if I have two or more categorical outcome variables? What type of analyses are available to me?
For example If I have data for School Grade [A, B, C, D, F], HS Drop Out[yes,no] Unemployed [yes,no], and I want to predict these three outcomes with Kindergarten Reading Ach. [0-100], Kindergarten Math Ach. [0-100], and Family SES [0-100]. Ignoring the validity of my model. How would I analyze those three categorical outcomes variables with my three continuous predictor variables?
 A: You are unclear on whether you are interested in the correlation structure of the outcome variables, or if you just want to handle polytomous data in a coherent way. I will assume the latter since it is more fundamental to the problem.
Since you have an ordered categorical outcome variable (grade), the easiest and most intuitive solution is probably proportional-odds logistic regression. This avoids having to impose an explicit continuous mapping (which often involves making an assumption, e.g. linearity, or doing a hack). The concept is poorly explained on the wikipedia page http://en.wikipedia.org/wiki/Ordered_logit and further reading is available in Agresti's Categorical Data Analysis (note, it is NOT in Introduction to Categorical Data Analysis which is basically a slimmed-down version of the former), or on http://www.stat.uchicago.edu/~pmcc/reports/prop_odds.pdf
It is implemented in R as the function polr() in the library MASS. You can type 'library(MASS)' followed by '?polr' for help.
A quick explanation would not go amiss. Note that logistic regression is basically defined by the model $\mathrm{logit}(pr(Y=1\,\big|\,X=x))=\beta_0+\sum_i^p \beta_i x_i$ under the usual assumptions of conditional independence.
Instead, let's take this model and say $$\mathrm{logit}(pr(\mathrm{Grade}>F\,\big|\,X=x))=\beta_{0F}+\sum_i^p \beta_i x_i.$$ This is equivalent to defining $Y=\chi(\mathrm{Grade}>F)$ and doing logistic regression. Now we'll continue and consider the conditional probability $$\mathrm{logit}(pr(\mathrm{Grade}>D\,\big|\,\mathrm{Grade}>F,X=x))=\beta_{0D}+\sum_i^p \beta_i x_i.$$ This means that the contribution of the explanatory variables (i.e. $\beta_i$) stays the same, but since $\big|\beta_{0D}\big|<\infty$, getting a C or better is harder than getting a D or better, in the sense that for any value of $x$, the probability is lower.
Obviously, repeating this for $pr(\mathrm{Grade}>C\,\big|\,\mathrm{Grade}>D, X=x)$ and so on will specify the model completely. In the end, you have four intercepts ($\beta_{0F},\beta_{0D},\beta_{0C},\beta_{0B}$) instead of just one ($\beta_0$).
Intuitively, a student's "quality", plus a random unobserved factor, must be high enough to overcome the hurdle of getting a C. The proportional odds logistic regression model basically says that it's probability $\mathrm{invlogit}(-\beta_{0D})$ to get from a $D$ to a $C$, no matter what quality $x$ one has. This is the assumption. Details are in the University of Chicago link above.
A: The idea of having different estimators for each output, as suggested by the other answers is part of the answer.  My suggestion would be to add more detail to the modeling of the system -- particularly in terms of the conditional independence relationships.  
For the groups that are conditionally independent, simply build separate models for each grup.  For the groups that are dependent, define a variables -- one for each group -- that live in the cartesian product space of the variables in question.  Estimates of this single variable can be mapped back to the multiple categorical variables in question.  Then everything works using the standard tool for each group.
In your case, the full product space of Grade x Drop_out x Unemployed would have 20 categories.
A: Nobody's submitted a direct response to this question yet IMO:
Timothy basically just suggested a univariate ordered logistic regression, which doesn't address the multivariate nature of the model you're interested in (not helpful).
Sifta suggests breaking the multivariate problem down into univariate analyses (not helpful).
Sharshofski starts talking about you converting previously continuous variables into categorical ones (not helpful).
McGuire seems unfamiliar with MANCOVA (not helpful).
So I'll share what I found online:
In principle, three categorical variables with 5, 2 and 2 levels will define a single categorical variable with 5 x 2 x 2 = 20 levels, so in principle, you can reduce the problem to one with a single categorical response and consider all sorts of techniques, such as multinomial models, for example.  There are heaps more of them, of course.  No one would ever go that way, though.  In practice you would separately predict each component using standard univariate models.  This is precisely how you fit ordinary multivariate linear regression models, of course.
This is an old problem and the classic text is "Discrete Multivariate Analysis" by Bishop, Fienberg and Holland (1st ed. in 1974!).  The (slightly) more recent book by Fahrmeir and Tutz on "Multivariate Statistical Analysis based on Generalized Linear Models" (Springer, 2001) is probably a bit more accessible, though be warned, however you look at it, this is tough territory.﻿
A: What is the format of your output before you convert it to a category? Since your inputs are continuous, I assume your model produces a continuous output at some point before you map that to the ultimate category.
When there are only two possible categories and your model is essentially a classifier, you can take a look at the ROC curve for your model. That'll allow you to pick a classifying threshold which maximizes the characteristics you want, i.e. minimizing false positives, false negatives, or the number of misclassified points in your data.
In my experience when I've had to deal with multi-class (multi-category in your words) modeling, I've trained a separate classifier for each class, then combined all the classifiers by choosing the one with the highest confidence to determine the class of a particular data point. See the one-vs-rest algorithm here. You should still be able to tune the individual classifiers using the ROC curve
A: Multiple outcomes (Y1, Y2, Y3) disqualifies regression. Structural equation modeling (SEM) is suited for multiple outcomes. However, if prediction of multiple dependent variables is needed, then you probably want to use causality, or more specifically, the structural approach to causation. Counterfactuals and logical approaches to causation and causal relationships will aid your model design, specification, confirmation and analysis also. 
See Judea Pearl (2009) Causality: Models, reasoning, and inference; 
Stephen L. Morgan and Christopher Winship (2007) Counterfactuals and causal inference: Methods and principles for social research; and 
Bill Shipley (2000) Cause and Correlation in Biology: A user's guide to Path Analysis, Structural Equations and Causal Inference. 
Best wishes.
