Simulation de régression linéaire multiple

14

Je suis nouveau dans le langage R. Je voudrais savoir comment simuler à partir d'un modèle de régression linéaire multiple qui remplit les quatre hypothèses de la régression.


D'accord, merci.

Disons que je veux simuler les données sur la base de cet ensemble de données:

y<-c(18.73,14.52,17.43,14.54,13.44,24.39,13.34,22.71,12.68,19.32,30.16,27.09,25.40,26.05,33.49,35.62,26.07,36.78,34.95,43.67)
x1<-c(610,950,720,840,980,530,680,540,890,730,670,770,880,1000,760,590,910,650,810,500)
x2<-c(1,1,3,2,1,1,3,3,2,2,1,3,3,2,2,2,3,3,1,2)

fit<-lm(y~x1+x2)
summary(fit)

alors j'obtiens la sortie:

Call:
lm(formula = y ~ x1 + x2)

Residuals:
     Min       1Q   Median       3Q      Max 
-13.2805  -7.5169  -0.9231   7.2556  12.8209 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)   
(Intercept) 42.85352   11.33229   3.782  0.00149 **
x1          -0.02534    0.01293  -1.960  0.06662 . 
x2           0.33188    2.41657   0.137  0.89238   
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1

Residual standard error: 8.679 on 17 degrees of freedom
Multiple R-squared:  0.1869,    Adjusted R-squared:  0.09127 
F-statistic: 1.954 on 2 and 17 DF,  p-value: 0.1722

Ma question est de savoir comment simuler de nouvelles données qui imitent les données originales ci-dessus?

Ni Hisham Haron
la source

Réponses:

28
  1. Si vous ne les avez pas déjà, commencez par configurer des prédicteurs, , , ...x 2x1x2

  2. Choisissez les coefficients de population («vrais») de vos prédicteurs, les , y compris , l'ordonnée à l'origine.β 0βiβ0

  3. Choisissez la variance d'erreur, ou de manière équivalente sa racine carrée, σσ2σ

  4. générer le terme d'erreur, , comme vecteur normal aléatoire indépendant, avec une moyenne de 0 et une variance σ 2εσ2

  5. Soity=β0+β1x1+β2x2+...+βkxk+ε

alors vous pouvez régressent l' sur vos s »xyx

Par exemple, dans R, vous pourriez faire quelque chose comme:

x1 <- 11:30
x2 <- runif(20,5,95)
x3 <- rbinom(20,1,.5)

b0 <- 17
b1 <- 0.5
b2 <- 0.037
b3 <- -5.2
sigma <- 1.4

eps <- rnorm(x1,0,sigma)
y <- b0 + b1*x1  + b2*x2  + b3*x3 + eps

produit une seule simulation de partir du modèle. Puis en cours d'exécutiony

 summary(lm(y~x1+x2+x3))

donne

Call:
lm(formula = y ~ x1 + x2 + x3)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.6967 -0.4970  0.1152  0.7536  1.6511 

Coefficients:
            Estimate Std. Error t value Pr(>|t|)    
(Intercept) 16.28141    1.32102  12.325 1.40e-09 ***
x1           0.55939    0.04850  11.533 3.65e-09 ***
x2           0.01715    0.01578   1.087    0.293    
x3          -4.91783    0.66547  -7.390 1.53e-06 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1   1

Residual standard error: 1.241 on 16 degrees of freedom
Multiple R-squared:  0.9343,    Adjusted R-squared:  0.9219 
F-statistic: 75.79 on 3 and 16 DF,  p-value: 1.131e-09

Vous pouvez simplifier cette procédure de plusieurs façons, mais je pensais que le préciser aiderait au début.

Si vous souhaitez simuler un nouveau aléatoire mais avec les mêmes coefficients de population, il suffit de réexécuter les deux dernières lignes de la procédure ci-dessus (générer un nouveau aléatoire et ), correspondant aux étapes 3 et 4 de l'algorithme.yepsy

Glen_b -Reinstate Monica
la source
Est-il possible de modifier l'erreur-type des estimations? J'ai utilisé un script légèrement modifié ( rnorm()au lieu de 11:30), mais peu importe à quel point j'augmente l'erreur (sigma), les erreurs standard de l'estimation sont à peu près similaires.
Daniel
2

Voici un autre code pour générer une régression linéaire multiple avec des erreurs suivant la distribution normale:

sim.regression<-function(n.obs=10,coefficients=runif(10,-5,5),s.deviation=.1){

  n.var=length(coefficients)  
  M=matrix(0,ncol=n.var,nrow=n.obs)

  beta=as.matrix(coefficients)

  for (i in 1:n.var){
    M[,i]=rnorm(n.obs,0,1)
  }

  y=M %*% beta + rnorm(n.obs,0,s.deviation)

  return (list(x=M,y=y,coeff=coefficients))

}
TPArrow
la source