source("http://192.38.117.59/~linearpredictors/datafiles/readVitaminD.R")

women <- subset(vitaminD, category == 2)
## Poland is now the reference level
women$country <- relevel(women$country, ref = "6")
women$bmigroup <- cut(women$bmi, c(0, 25, 30, Inf),
                            right = F,
                            labels = c("Normal", "Slight","Obese")) 

countryonly <- lm(log10(vitd) ~ country, data = women)
bmionly <- lm(log10(vitd) ~ bmigroup, data = women)
adjusted <- lm(log10(vitd) ~ bmigroup + country, data = women)