Influence of two factors on non-normal outcome. Which test?

I measured two variables (outcomes) $Y_a$ and $Y_b$ for several subjects. $Y_a$ and $Y_b$ are continuous variables defined on $[0,100]$ and have a bimodal distribution. I tried using a transformation but came up with nothing worthy. I have two independent variables $X_1$ (factor with 2 levels) and $X_2$ (factor with 3 levels) characterizing each subject. I would like to investigate if $X_1$ and $X_2$ influence $Y_a$ or $Y_b$. In normally distributed outcomes I would have used lm(Ya ~ X1 + X2). What is the non-parametric equivalent test? Is it the Kruskal-Wallis? Do I have to recode $X_a$ and $X_b$ in one variable?

Linear models do not make assumptions about the distribution of the dependent variable, they make assumptions about the distribution of the error, as measured by the residuals. That said, if the Y variables are bimodal, you may want to think about quantile regression. I wrote about how to do this in SAS but that paper also shows some basic points not specific to SAS. Since evidently you use R you can look into the quantreg package.

I suggest quantile regression because, if Y is bimodal, it seems likely to me that different things may be related to changes in Y at or near the different modes.

I am not sure what your last sentence refers to.

@Peter Flom's point on precisely what is assumed is crucial.

In addition, much depends on the nature of the bimodality and how marked it is. You don't say what transformations you tried, but the most usual transformations such as logarithm or square root are indeed unlikely to help here. However, If the modes are near 0 and near 100 so that your distribution is approximately U-shaped, then a logit transformation should help, once you divide by 100, or equivalently generalise the transformation to match the bounds.

For example, a logit transformation always removes bimodality from a beta distribution.

However, you may have exact zeros or exact unities, and for that reason and on other grounds, a generalized linear model with logit link is the better way to go. See also http://www.stata-journal.com/sjpdf.html?articlenum=st0147 That paper is geared to Stata implementation, but the principles are as usual more general.

Note that the Kruskal-Wallis test's way of coping with bimodality is to ignore it. Only in that sense is it a solution. Furthermore, it cannot handle two predictors or interaction. Recoding two predictors to a composite predictor, which seems to be what you are asking about, again solves that problem only by ignoring it.