I fit a model to my data with the following formula after a stepwise selection (x1 to x5 stand for my variables):
lm(formula = outcome ~ x1+ x3 + x4 + x2+ x5+ poly(x4, 2) + poly(x2, 2) + poly(x5, 2) + poly(x1, 2) + x1:x2 + x1:x4+ x1:x3 + x4:x2 + x4:x5 + x3:x4+ x3:x2 + x2:x5 + x3:x5, data = data)
VIF shows strong multicollinearity in the interaction terms and one/two poly terms. Removing predictors with very high VIF increases the model error. I have heard that varpart() function in R can handle the multicollinearity issue, but that can only accept maximum 4 terms, which is below the number of high VIF terms in my model. Is multicollinearity an issue in this case and might be due to its polynomial nature and can be ignored or how best to handle it? The variables are all related biologically, but the VIF values are around 1 in a simple linear model of factors.
Thanks for any help