#Simulation
set.seed(1560101)
n<-200
p<-100
Sig<-rWishart(1,p+1,diag(1,p))[,,1]
library(MASS)
X<-mvrnorm(n,mu=rep(0,p),Sigma=Sig)

beta<-rep(0,p+1)
beta[1]<-10
beta[2:6]<-c(-0.3,-2,1,-2,0.5)
Xmat<-cbind(rep(1,n),X)

sig<-5.0

mean.Y<-Xmat %*% beta

Y<-rnorm(n)*sig+mean.Y

index<-1:p

############################################################
#No penalty
alpha<-0.05/p
selected.ls<-index[coef(summary(lm(Y~X)))[-1,4]<alpha]

library(MASS)

