Tirez sur les valeurs de p et r au carré d’une régression linéaire

Question

Comment extrayez-vous la valeur p (pour la signification du coefficient de la variable explicative unique non nul) et la valeur R au carré d'un modèle de régression linéaire simple? Par exemple...

 x = cumsum(c(0, runif(100, -1, +1)))
y = cumsum(c(0, runif(100, -1, +1)))
fit = lm(y ~ x)
summary(fit)
 

Je sais que summary(fit) affiche la valeur p et la valeur R au carré, mais je veux pouvoir les coller dans d'autres variables.

Answer 1

Vous pouvez retourner le r-carré de la valeur directement à partir du résumé d'objet summary(fit)$r.squared. Voir names(summary(fit)) pour une liste de tous les articles que vous pouvez extraire directement.

Ce blog présente une fonction pour renvoyer la valeur de p:

lmp <- function (modelobject) {
    if (class(modelobject) != "lm") stop("Not an object of class 'lm' ")
    f <- summary(modelobject)$fstatistic
    p <- pf(f[1],f[2],f[3],lower.tail=F)
    attributes(p) <- NULL
    return(p)
}

> lmp(fit)
[1] 1.622665e-05

Sinon, vous pouvez saisir la valeur de p à partir de l' anova(fit) objet dans un mode similaire à l'objet de résumé ci-dessus.

Answer 2

Notez que summary(fit) génère un objet avec toutes les informations dont vous avez besoin. Les vecteurs bêta, se, t et p y sont stockés. Obtenez les valeurs p en sélectionnant la 4ème colonne de la matrice de coefficients (stockée dans l'objet récapitulatif):

 summary(fit)$coefficients[,4] 
summary(fit)$r.squared
 

Essayez str(summary(fit)) pour voir toutes les informations que cet objet contient.

Edit: J'avais mal interprété la réponse de Chase qui vous dit en gros comment arriver à ce que je donne ici.

Answer 3

Vous pouvez voir la structure de l'objet renvoyé par summary() en appelant str(summary(fit)) . Chaque pièce est accessible avec $ . La valeur p de la statistique F provient plus facilement de l'objet renvoyé de anova .

En bref, vous pouvez faire ceci:

 rSquared <- summary(fit)$r.squared
pVal <- anova(fit)$'Pr(>F)'[1]
 

Answer 4

Alors que les deux réponses ci-dessus sont bonnes, la procédure pour extraire des parties d'objets, d'ordre plus général.

Dans de nombreux cas, les fonctions renvoient des listes, et les composants individuels peuvent être accédées en utilisant str() ce qui permet d'imprimer les composants ainsi que leurs noms. Vous pouvez ensuite y accéder à l'aide de l'opérateur, c'est à dire myobject$componentname.

Dans le cas de lm objets, il ya un certain nombre de méthodes prédéfinies que l'on peut utiliser comme coef(), resid(), summary() etc, mais vous ne pourrez pas toujours être de la chance.

Answer 5

J'ai utilisé ce lmp fonction tout à fait un grand nombre de fois.

Et à un moment, j'ai décidé d'ajouter de nouvelles fonctionnalités pour améliorer l'analyse des données. Je ne suis pas un expert dans R ou des statistiques, mais les gens sont généralement à la recherche à des informations différentes d'une régression linéaire :

• la p-valeur
• a et b
• r2
• et bien sûr, l'aspect de la distribution des points

Prenons un exemple. Vous avez ici

Voici un exemple reproductible avec des variables différentes:

Ex<-structure(list(X1 = c(-36.8598, -37.1726, -36.4343, -36.8644, 
-37.0599, -34.8818, -31.9907, -37.8304, -34.3367, -31.2984, -33.5731
), X2 = c(64.26, 63.085, 66.36, 61.08, 61.57, 65.04, 72.69, 63.83, 
67.555, 76.06, 68.61), Y1 = c(493.81544, 493.81544, 494.54173, 
494.61364, 494.61381, 494.38717, 494.64122, 493.73265, 494.04246, 
494.92989, 494.98384), Y2 = c(489.704166, 489.704166, 490.710962, 
490.653212, 490.710612, 489.822928, 488.160904, 489.747776, 490.600579, 
488.946738, 490.398958), Y3 = c(-19L, -19L, -19L, -23L, -30L, 
-43L, -43L, -2L, -58L, -47L, -61L)), .Names = c("X1", "X2", "Y1", 
"Y2", "Y3"), row.names = c(NA, 11L), class = "data.frame")


library(reshape2)
library(ggplot2)
Ex2<-melt(Ex,id=c("X1","X2"))
colnames(Ex2)[3:4]<-c("Y","Yvalue")
Ex3<-melt(Ex2,id=c("Y","Yvalue"))
colnames(Ex3)[3:4]<-c("X","Xvalue")

ggplot(Ex3,aes(Xvalue,Yvalue))+
          geom_smooth(method="lm",alpha=0.2,size=1,color="grey")+
          geom_point(size=2)+
          facet_grid(Y~X,scales='free')


#Use the lmp function

lmp <- function (modelobject) {
  if (class(modelobject) != "lm") stop("Not an object of class 'lm' ")
  f <- summary(modelobject)$fstatistic
    p <- pf(f[1],f[2],f[3],lower.tail=F)
    attributes(p) <- NULL
    return(p)
    }

# create function to extract different informations from lm

lmtable<-function (var1,var2,data,signi=NULL){
  #var1= y data : colnames of data as.character, so "Y1" or c("Y1","Y2") for example
  #var2= x data : colnames of data as.character, so "X1" or c("X1","X2") for example
  #data= data in dataframe, variables in columns
  # if signi TRUE, round p-value with 2 digits and add *** if <0.001, ** if < 0.01, * if < 0.05.

  if (class(data) != "data.frame") stop("Not an object of class 'data.frame' ")
  Tabtemp<-data.frame(matrix(NA,ncol=6,nrow=length(var1)*length(var2)))
  for (i in 1:length(var2))
       {
  Tabtemp[((length(var1)*i)-(length(var1)-1)):(length(var1)*i),1]<-var1
  Tabtemp[((length(var1)*i)-(length(var1)-1)):(length(var1)*i),2]<-var2[i]
  colnames(Tabtemp)<-c("Var.y","Var.x","p-value","a","b","r^2")

  for (n in 1:length(var1))
  {
  Tabtemp[(((length(var1)*i)-(length(var1)-1))+n-1),3]<-lmp(lm(data[,var1[n]]~data[,var2[i]],data))

  Tabtemp[(((length(var1)*i)-(length(var1)-1))+n-1),4]<-coef(lm(data[,var1[n]]~data[,var2[i]],data))[1]

  Tabtemp[(((length(var1)*i)-(length(var1)-1))+n-1),5]<-coef(lm(data[,var1[n]]~data[,var2[i]],data))[2]

  Tabtemp[(((length(var1)*i)-(length(var1)-1))+n-1),6]<-summary(lm(data[,var1[n]]~data[,var2[i]],data))$r.squared
  }
  }

  signi2<-data.frame(matrix(NA,ncol=3,nrow=nrow(Tabtemp)))
  signi2[,1]<-ifelse(Tabtemp[,3]<0.001,paste0("***"),ifelse(Tabtemp[,3]<0.01,paste0("**"),ifelse(Tabtemp[,3]<0.05,paste0("*"),paste0(""))))
  signi2[,2]<-round(Tabtemp[,3],2)
  signi2[,3]<-paste0(format(signi2[,2],digits=2),signi2[,1])

  for (l in 1:nrow(Tabtemp))
    {
  Tabtemp$"p-value"[l]<-ifelse(is.null(signi),
         Tabtemp$"p-value"[l],
         ifelse(isTRUE(signi),
                paste0(signi2[,3][l]),
                Tabtemp$"p-value"[l]))
  }

   Tabtemp
}

# ------- EXAMPLES ------

lmtable("Y1","X1",Ex)
lmtable(c("Y1","Y2","Y3"),c("X1","X2"),Ex)
lmtable(c("Y1","Y2","Y3"),c("X1","X2"),Ex,signi=TRUE)

Il y a certainement une solution plus rapide que de cette fonction, mais il fonctionne.