Lineare Regression und Gruppierung nach in R

Question

Ich möchte eine lineare Regression in R mit der lm()-Funktion durchführen. Meine Daten sind eine jährliche Zeitreihe mit einem Feld für das Jahr (22 Jahre) und einem anderen für den Bundesstaat (50 Bundesstaaten). Ich möchte eine Regression für jeden Bundesstaat anpassen, sodass ich am Ende einen Vektor von lm-Antworten habe. Ich kann mir vorstellen, für jede Zustandsschleife eine Schleife zu machen, dann die Regression innerhalb der Schleife durchzuführen und die Ergebnisse jeder Regression zu einem Vektor hinzuzufügen. Das scheint jedoch nicht sehr R-typisch zu sein. In SAS würde ich eine 'by'-Anweisung verwenden und in SQL würde ich ein 'group by' machen. Wie funktioniert das in R?

Answer 1

Ich denke, es lohnt sich, den purrr::map Ansatz zu diesem Problem hinzuzufügen.

library(tidyverse)

d <- data.frame(state=rep(c('NY', 'CA'), c(10, 10)),
                                 year=rep(1:10, 2),
                                 response=c(rnorm(10), rnorm(10)))

d %>% 
  group_by(state) %>% 
  nest() %>% 
  mutate(model = map(data, ~lm(response ~ year, data = .)))

Siehe die Antwort von @Paul Hiemstra für weitere Ideen zur Verwendung des broom Pakets mit diesen Ergebnissen.

Answer 2

## Erstellen von Fake-Daten
 ngroups <- 2
 group <- 1:ngroups
 nobs <- 100
 dta <- data.frame(group=rep(group,each=nobs),y=rnorm(nobs*ngroups),x=runif(nobs*ngroups))
 head(dta)
#--------------------
  group          y         x
1     1  0.6482007 0.5429575
2     1 -0.4637118 0.7052843
3     1 -0.5129840 0.7312955
4     1 -0.6612649 0.9028034
5     1 -0.5197448 0.1661308
6     1  0.4240346 0.8944253
#------------ 
## Funktion zur Extraktion der Ergebnisse eines Modells
 foo <- function(z) {
   ## Koeffizienten und Standardfehler in einem Datenrahmen
   mr <- data.frame(coef(summary(lm(y~x,data=z))))
   ## Fügen Sie Zeilenbezeichnungen hinzu (Prädiktoren/Unabhängige Variablen)
   mr$predictor <- rownames(mr)
   mr
 }
 ## Schauen Sie, ob es funktioniert
 foo(subset(dta,group==1))
#=========
              Estimate Std..Error   t.value  Pr...t..   predictor
(Intercept)  0.2176477  0.1919140  1.134090 0.2595235 (Intercept)
x           -0.3669890  0.3321875 -1.104765 0.2719666           x
#----------
## Eine Option: Verwenden Sie den Befehl by
 res <- by(dta,dta$group,foo)
 res
#=========
dta$group: 1
              Estimate Std..Error   t.value  Pr...t..   predictor
(Intercept)  0.2176477  0.1919140  1.134090 0.2595235 (Intercept)
x           -0.3669890  0.3321875 -1.104765 0.2719666           x
------------------------------------------------------------ 
dta$group: 2
               Estimate Std..Error    t.value  Pr...t..   predictor
(Intercept) -0.04039422  0.1682335 -0.2401081 0.8107480 (Intercept)
x            0.06286456  0.3020321  0.2081387 0.8355526           x

## Die Verwendung des Pakets plyr ist besser
 library(plyr)
 res <- ddply(dta,"group",foo)
 res
#----------
  group    Estimate Std..Error    t.value  Pr...t..   predictor
1     1  0.21764767  0.1919140  1.1340897 0.2595235 (Intercept)
2     1 -0.36698898  0.3321875 -1.1047647 0.2719666           x
3     2 -0.04039422  0.1682335 -0.2401081 0.8107480 (Intercept)
4     2  0.06286456  0.3020321  0.2081387 0.8355526           x

Answer 3

Ich weiß, meine Antwort kommt ein wenig spät, aber ich habe nach einer ähnlichen Funktionalität gesucht. Es scheint, dass die integrierte Funktion 'by' in R auch die Gruppierung einfach durchführen kann:

?by enthält das folgende Beispiel, das pro Gruppe passt und die Koeffizienten mit sapply extrahiert:

require(stats)
## angenommen, wir möchten die Koeffizienten nach Gruppe extrahieren
tmp <- with(warpbreaks,
            by(warpbreaks, tension,
               function(x) lm(breaks ~ wool, data = x)))
sapply(tmp, coef)

Answer 4

Die lm() Funktion oben ist ein einfaches Beispiel. Übrigens, ich stelle mir vor, dass Ihre Datenbank die Spalten in folgender Form hat:

year state var1 var2 y...

Meiner Meinung nach können Sie den folgenden Code verwenden:

require(base) 
library(base) 
attach(data) # data = Ihre Datenbank
             #state ist Ihr Label für die Spalten der Bundesländer
modell<-by(data, data$state, function(data) lm(y~I(1/var1)+I(1/var2)))
summary(modell)

Answer 5

Die Frage scheint zu sein, wie Regressionsfunktionen mit Formeln aufgerufen werden können, die innerhalb einer Schleife modifiziert werden.

So können Sie es tun (unter Verwendung des Diamantdatensatzes):

attach(ggplot2::diamonds)
strCols = names(ggplot2::diamonds)

formel <- list(); modell <- list()
for (i in 1:1) {
  formel[[i]] = paste0(strCols[7], " ~ ", strCols[7+i])
  modell[[i]] = glm(formel[[i]]) 

  # dann können Sie die Ergebnisse plotten oder etwas anderes tun ...
  png(filename = sprintf("diamonds_price=glm(%s).png", strCols[7+i]))
  par(mfrow = c(2, 2))      
  plot(modell[[i]])
  dev.off()
  }