Are there any circumstances where stepwise regression should be used? Stepwise regression had been overused in many biomedical papers in the past but this appears to be improving with better education of its many issues. Many older reviewers however do still ask for it. What are the circumstances where stepwise regression has a role and should be used, if any?
 A: I am not aware of situations, in which stepwise regression would be the preferred approach. It may be okay (particularly in its step-down version starting from the full model) with bootstrapping of the whole stepwise process on extremely large datasets with $n>>p$. Here $n$ is the number of observations in an continuous outcome (or number of records with an event in survival analysis) $p$ is the number of candidate predictors including all considered interactions - i.e. when any even small effects become very clear and it does not matter so much how your do your model building (that would mean that $n$ would be much larger than $p$ than by substantially more than the sometimes quoted factor of 20).
Of course the reason most people are tempted to do something like stepwise regression is, 


*

*because it is not computationally intensive (if you do not do the proper bootstrapping, but then your results are pretty unreliable), 

*because it provides clear cut "is in the model" versus "is not in the model" statements (which are very unreliable in standard stepwise regression; something that proper bootstrapping will usually make clear so that these statements will usually not be so clear) and 

*because often $n$ is smaller, close to or just a bit larger than $p$. 


I.e. a method like stepwise regression would (if it had good operating characteristics) be especially attractive in those situations, when it does not have good operating characteristics.
A: Two cases in which I would not object to seeing step-wise regression are


*

*Exploratory data analysis

*Predictive models


In both these very important use cases, you are not so concerned about traditional statistical inference, so the fact that p-values, etc., are no longer valid is of little concern.
For example, if a research paper said "In our pilot study, we used step-wise regression to find 3 interesting variables out of 1000. In a follow-up study with new data, we showed these 3 interesting variables were strongly correlated with the outcome of interest", I would have no problem with the use of step-wise regression. Similarly, "We used step-wise regression to build a predictive model. This out-preformed alternative model X in our hold-out data set in regards to MSE" is totally fine with me as well. 
To be clear, I am not saying that step-wise regression is the best way to approach these problems. But it is easy and may give you satisfactory solutions. 
EDIT:
In the comments, there is a question of whether stepwise AIC can actually be useful for prediction. Here's a simulation that shows it doing much better than linear regression with all the covariates, and nearly as well as elastic nets with the penalty chosen by cross-validation. 
I wouldn't take this simulation as the end of the discussion; it's not too hard to come up with a scenario in which step-wise AIC will preform worse. But it's really not an unreasonable scenario, and exactly the type of situation that elastic nets are designed for (high correlation of covariates with very few large effects)!  
library(leaps)
library(glmnet)
nRows <- 1000
nCols <- 500

# Seed set For reproducibility. 
# Try changing for investigation of reliability of results
set.seed(1)

# Creating heavily correlated covariates
x_firstHalf  <- matrix(rnorm(nRows * nCols / 2), nrow = nRows)
x_secondHalf <- x_firstHalf + 0.5 * 
                matrix(rnorm(nRows * nCols / 2), nrow = nRows) 
x_mat        <- cbind(x_firstHalf, x_secondHalf) + rnorm(nRows)

# Creating beta's. Most will be of very small magnitude
p_large = 0.01
betas <- rnorm(nCols, sd = 0.01) + 
         rnorm(nCols, sd = 4) * rbinom(nCols, size = 1, prob = p_large)
y     <- x_mat %*% betas + rnorm(nRows, sd = 4)

all_data           <- data.frame(y, x_mat)
colnames(all_data) <- c('y', paste('x', 1:nCols, sep = '_'))

# Holding out 25% of data for validation
holdout_index <- 1:(nRows * .25) 
train_data    <- all_data[-holdout_index, ]
validate_data <- all_data[holdout_index, ]

mean_fit <- lm(y ~ 0, data = train_data)
full_fit <- lm(y ~ ., data = train_data)
step_fit <- step(mean_fit, 
                 scope = list(lower = mean_fit, upper = full_fit), 
                 direction = "forward", steps = 20, trace = 0)

glmnet_cvRes <- cv.glmnet(x = as.matrix(train_data[,-1]), 
                          y = as.numeric(train_data$y)   )

full_pred   <- predict(full_fit, validate_data)
step_pred   <- predict(step_fit, validate_data)
glmnet_pred <- predict(glmnet_cvRes, as.matrix(validate_data[,-1]), s='lambda.min')

sd(full_pred - validate_data$y)    # [1] 6.426117
sd(step_pred - validate_data$y)    # [1] 4.233672
sd(glmnet_pred - validate_data$y)  # [1] 4.127171
# Note that stepwise AIC does considerably better than using all covariates 
# in linear regression, and not that much worse than penalized methods
# with cross validation!!

Side note:
I'm really not a fan of stepwise regression for many, many reasons, so I feel somewhat awkward having taken this stance in defense of it. But I merely think it's important to be precise about exactly what I don't like about it. 
