Combined antiretroviral therapies (CART) can successfully suppress HIV in the serum and bring its viral load below detection rate. However, drug resistance remains a major challenge. As resistance patterns vary between patients, personalized therapy is required. Automatic systems for therapy personalization exist and were shown to better predict therapy outcome than HIV experts in some settings. However, these systems focus only on selecting the therapy most likely to suppress the virus for several weeks, a choice that may be suboptimal over the longer term due to evolution of drug resistance. We present a novel generative model for HIV drug resistance evolution. This model is based on factorial HMMs, applying a novel collapsed Gibbs Sampling algorithm for approximate learning. Using the suggested model, we obtain better therapy outcome predictions than existing methods and recommend therapies that may be more effective over the long term. We demonstrate our results using simulated data and using real data from the EuResist dataset.