basic procedure
Regression, logistic
Generalized linear model (example is a case-control study)
mod = glm(response_var~covariate_va1 + covariate_var2 + interaction_var1:interaction_var2, data = dat, family = binomial)
summary(mod) # gives coefficients and p values
exp(mod$coefficients) # exponentiate to get odds ratio
exp(confint(mod)) # exponentiate the confidence interval
Regression, variable selection : https://www.r-bloggers.com/variable-selection-using-automatic-methods/