# R - Identifying the most statistically significant variables that determine a target

My dataframe looks something like this -

df <- data.frame(Speed_Satisfaction=sample(seq(from = 0, to = 1, by = 0.01), size = 1000, replace = TRUE),
Cleanliness_Satisfaction=sample(seq(from = 0, to = 1, by = 0.01), size = 1000, replace = TRUE),
Efficiency_Satisfaction=sample(seq(from = 0, to = 1, by = 0.01), size = 1000, replace = TRUE),
Variety_Satisfaction=sample(seq(from = 0, to = 1, by = 0.01), size = 1000, replace = TRUE),
Friendliness_Satisfaction=sample(seq(from = 0, to = 1, by = 0.01), size = 1000, replace = TRUE),
Overall_Satisfaction=sample(seq(from = 0, to = 1, by = 0.01), size = 1000, replace = TRUE))


My Objective

I'd like to determine which of my 5 variables is most important in determining Overall_Satisfaction for services provided in a restaurant. And rank these variables in order of their importance.

The satisfaction scores in my df are on a scale from 0 to 1, where 1 equals 100% satisfaction and 0 is 0% satisfaction.

My Approach

I chose to model my data using multiple linear regression (the lm function in R) to identify the most significant variables.

fit <- lm(Overall_Satisfaction ~ Speed_Satisfaction + Cleanliness_Satisfaction + Variety_Satisfaction + Friendliness_Satisfaction, data=df)


And the analyse the results using summary

This way, I can study the Estimates (negative or positive indicating the direction of the relationship with the target variable) and the P value to rank variables according to how statistically significant they are.

My Question

Is this the best approach for the problem I'm trying to solve with the data I have, or are there other methods that yield better results?

I am still exploring statistics and R, so your inputs will be greatly appreciated.

As John Coleman said in his comment, Cross Validated might be a better place to look for an answer, and you should ask there for more details. Also, in your idealized example above, most problems below don't apply since all the variables seem to have the same distribution, however you need to keep them in mind. For now, the main points would be:

1) In general, don't look at the estimates, because they are a function of the units in which they are measured. So for example, if you measure something in kilos and in grams, the latter coefficient will be 1000 larger, but without it being a better predictor. That being said, your variables are on the same scale, so the problem goes away.

2) p-values measure how likely you are to observe the estimated effect (or a more extreme one) under the null hypothesis. They are a function of your coefficient, the precision of the estimation (standard error), the sample size, and the null hypothesis. So again, not a great way to gauge importance.

Suggestions:

1) Look at standardized coefficients. In this case, you normalize your variables (subtract mean, divide by standard deviation) and the resulting coefficients are comparable. However, the interpretation changes. You can't say "a 1-unit increase in X, changes Y by ", but you can say "a 1-sd increase in X...". (In your example, the data are already standardized). You can implement this in R using the lm.beta package.

2) Remove the variables one by one, and see how the Adjusted changes. So if you remove , and the Adjusted decreases by , where , then you can say that is the most important predictor. However, the results are not necessarily reliable, and that's why there is an extensive literature on variable selection methods (backwards, forward, best subset, LASSOs, Ridge Regressions, Elastic Nets etc.). However, with the basic example you gave, a basic application of backwards selection should be fine.

This is not by any means a complete discussion, but I think it hits some very basic points you should keep in mind.

• While suggestion (2) is good to start with, I'd advise some caution in that your result will also depend on the order in which you remove variables. Incidentally, this is one of the reasons why stepwise regression methods can be very unreliable. May 2, 2018 at 20:55
• @MarcusCampbell, that's an excellent point! I was originally debating whether to get into it or not, and chose not to because of the simplicity of the OP's example. But I just edited the answer to capture the spirit of your comment. Thanks! May 2, 2018 at 21:09
• "p-values measure how likely you are to observe a coefficient under the null hypothesis" doesn't seem totally exact. I would say "p-valuez measure how likely you are to observe the at least the estimated effect under the null hypothesis" May 2, 2018 at 21:32
• @denis fair point! I edited the answer May 2, 2018 at 21:37

Variables such as satisfaction are usually best modeled with ordinal regression, e.g., the proportional odds model, because of the lack of interval scaling. But to your question, the data are ill-equipped for telling you the importance of variables in prediction, especially if the variables are correlated. I suggest using a method that is honest and exposes the difficulty of the task: using the bootstrap to compute 0.95 confidence limits for the ranks of the variables on their prediction ability. This can be based on the log-likelihood that is explained by the predictors (partial $\chi^2$ Wald tests; likelihood ratio tests better). A simple example of such a bootstrap exercise in a simulation of 12 predictors is in my RMS course notes.

If, when doing the bootstrap, one finds that the data are compatible with a certain variable being both the strongest and the weakest predictors, you have no basis for making a statement.

Background: look at a few regular bootstrap repetitions of model fits and see how the relative contributions of all the predictors are volatile.

I don't want to enter the whole discussion, I just want to point out a simple tool:

library(Hmisc)
rcorr(as.matrix(df))


It will give you the different correlations and the p-values associated. The strongest significant correlation would be a good start to determine your variable of interest