How can I deal with a covariate being defined only for a subset of my sample? My model looks like: ROI_size ~ diagnosis + medication_dose + sex + age. Specifically, I want to find the effect of disease (1 or 0), adjusted for current medication dose (measured in mg) on brain region size, while controlling for age and sex. Both of the latter variables are available for the entire dataset. However, the former (medication dose) is defined only for the patients (diagnosis = 1).
I have thought of two solutions. One would be to assign medication_dose = 0 to all healthy controls. However, this creates 1. collinearity with diagnosis 2. zero-inflated data and is per se not correct (I think), since the difference between 0 and 1 mg is not the same as the difference between 1 mg and 2 mg.
The other solution would be to adjust the ROI_size for medication_dose within the patient cohort, i.e., to run the model ROI_size ~ medication_dose and then by calculating its beta and then adding beta*medication_dose to the ROI_size for all patients.
My question in this case is if I should have age and sex in the second model, since including them will yield a more accurate beta for medication_dose. In that case, do I need to correct for age and sex separately for the healthy controls as well and then to find the effect of disease to just do ROI_size ~ Diagnosis?
Histogram for current medication use in patients

Histogram for current medication use in all



... there is no healthy control that has medication > 0 and no patient that has medication = 0.

That suggests a different two-step modeling approach from your second solution. Start without the medication.dose term in the model. That evaluates the critical issue of diagnosis directly, controlling for Sex and Age. Then add the medication.dose term and see whether it improves the fit; e.g., anova() comparison of those two nested models.
If the fit is not improved, then you have documented that adjusting for medication.dose doesn't matter (much) in terms of ROI_size when diagnosis, Sex and Age are taken into account.
If the fit is improved, then you have an estimate of the association between ROI_size and medication.dose. The model is then equivalent to what @whuber recommends in a different context for a predictor whose value is necessarily 0 for some cases. Here, the diagnosis predictor serves the role of the "loan indicator" in that question, with a corresponding interpretation of intercepts and regression coefficients.
There's an important caution about how to handle medication.dose, however. If you treat it as an untransformed continuous variable then you are imposing a linear association between medication.dose and ROI_size. If there are only a few different values of medication.dose then you might be best off treating it as a multi-level categorical predictor, perhaps ordinal. If there's a wide range of values, you might need to transform appropriately.
You would have a complete collinearity problem if there was only one value of medication.dose for all those with a diagnosis; in that case you would not be able to fit a model with both diagnosis and medication.dose uniquely. The partial collinearity that you fear with more than one medication.dose might tend to make the error estimates of the regression coefficients larger. If you are interested specifically in the association of medication.dose with ROI_size for those with diagnosis you might be better off restricting analysis to those with diagnosis. The question (which I don't think can be answered except by looking at the data) is whether having the correction for Age and Sex from extra cases without diagnosis makes up for that.
