Isolating a relationship's significance The thought experiment:
STEM = having a science, technology, engineering or mathematics degree (0/1), W = wages in the top quartile (0/1), LOC = individual location vector (integer values 1-20 representing different cities)

*

*When we regress W on STEM (reg w STEM), we find that the coefficient on STEM is positive

*When we regress W on STEM and LOC (reg W STEM LOC1 LOC2 ... LOC20), we find the coefficient on STEM is less positive

*When we regress STEM on LOC (reg STEM LOC1 LOC2 ... LOC20), we find that the coefficients on some LOC are positive and some are negative

*When we regress W and LOC (reg W LOC1 LOC2 ... LOC20), we find that the coefficients on some LOC are positive and some are negative

Q1) How can discuss the amount of the relationship between W = 1 and STEM = 1 that is due to STEM = 1 versus the fact that people with STEM = 1 may choose to live in locations (LOC) where W is more likely to be 1 and avoid locations where W is less likely to be 1?
Q2) What models, model output, and/or post estimation commands help create a more comprehensive story?
 A: Since your dependent variable (W) is Binary, you could consider the Logistic Regression approach. The Logistic Regression coefficients are in terms of log-odds and you can easily derive the odds or probabilities from these coefficients.
Logistic Regression allows you to estimate Marginal Effects. These can help you understand the effect of STEM on the probability of W (having wages in the top quartile). For example in Stata software, you can easily accomplish this using "mfx" command after logistic regression.
The interpretation of marginal effects is the change in probability of W if the individual has a STEM degree. Suppose the marginal effect of STEM is estimated to be 0.2, it means that the probability of having Wages in the top quartile increases by 0.2 if the individual has a STEM degree.
For locations, you could consider including dummy variables for locations. The upside of this is that you can estimate the marginal effects/probabilities and their effects on W for each location separately. The downside is that it will introduce a lot of coefficients in your model and restrict degrees of freedom. You will need a bigger dataset for this.
