Comparing different ways and advantages/disadvantages of summarizing/clustering and matching a set of explanatory variables? I have the following data. Feature 1... Feature 15 is a set of explanatory variables. I have a group variable Med_group to show the two different possible medication received by patients. And then a outcome variable Days_being_healthy as the number of days being healthy.
Feature1    Feature2    …   Feature 15      Med_group   Days_being_healthy
1           1               123             1            23
2           2               1234            1            534
4           3               345             1            65
1           2               12              2            85
5           2               3466            1            185
1           2               34              2            85
5           1               634             2            24

Since I believe the feature profile (#1-#15) is very different between the Med_group, I want to address for it. I am thinking there are three ways to do it:
1) Using t-test to compare Days_being_healthy between the two Med_group, while adjusting for all 15 features;
2) Similar to 1), but combine the 15 features into propensity score, and adjust the propensity score;
3) Doing an unsupervised clustering with the 15 features, then within each cluster, reveal the Med_group label, then compare the Days_being_healthy between the "1" and "2" groups.
I understand 2) has an advantage when there is a lot of features and few outcomes. But what might be the advantage of approach 3)?
 A: First, it would be a good idea to combine 1 and 2. Second, rather than including the propensity score in an adjusted t-test/regression, estimate balancing weights, and use those in a weighted t-test/regression that also adjusts for the features. Balancing weights are weights that when applied to your data set will make the weighted groups appear similar to each other on the features. A good method of weighting is generalized boosted modeling, but there are others.
In R, your command would look like: 
lm(Days_being_healthy ~ Med_group, data = data, weights = balancing.weights)

presuming your data is stored in data and your balancing weights are stored as a vector called balancing.weights. You can do a bootstrap of the entire process to estimate standard errors or a confidence interval. When choosing tuning parameters for boosting, I recommend "ks.max".
Your #3 would work but only on the condition that within each cluster, the groups look similar on the features. A common way to do this without machine learning is to estimate propensity scores and stratify in quantiles of the propensity score, then estimate treatment effects within each quantile and pool them. This approach has been shown to be more biased than weighting as I described above.
