Quiero realizar una de Regresión lineal utilizando los valores de p como criterio de selección, por ejemplo: en cada paso de caer variables que tienen el más alto es decir, el más insignificante de los valores de p, parando cuando todos los valores son importantes definidos por cierto umbral alfa.

Soy totalmente consciente de que debo usar el AIC (por ejemplo, el comando paso o stepAIC) o algún otro criterio en lugar de eso, pero mi jefe tiene ninguna comprensión de las estadísticas y de insistir en el uso de p-valores.

Si es necesario, puedo programar mi propia rutina, pero me pregunto si hay una que ya se implementó la versión de este.

  • Espero que tu jefe no lea stackoverflow. 😉
  • ¿Por qué exactamente de regresión paso a paso? ¿Cuántas variables se tienen?
  • Usted podría considerar la posibilidad de hacer esta pregunta aquí: stats.stackexchange.com
  • Yo también mi jefe no lea esto 😉 tengo 9 y variables a tener que probar de todo un poco, para dejar caer las variables «a mano» y el ajuste de un nuevo modelo mí es un poco mucho para escribir, así que me pregunto si hay una manera automatizada como con el «paso», sólo con p-valores.
  • Porque alguien ha utilizado la regresión paso a paso con los valores de p en el pasado (con el programa STATA supongo, pero no tenemos STATA más), y ella insiste en usar el mismo enfoque. Esa es la manera en que son los jefes … 😉
  • Podría ser más fácil de enseñar jefe de buenas estadísticas que a R a hacer malas estadísticas.
  • Sólo elige tres variables al azar – usted probablemente va a hacer igual de bien regresión escalonada.
  • ¿Su jefe también dígale a su médico qué medicamentos a prescribir y su mecánico de cómo arreglar su coche?

InformationsquelleAutor DainisZ | 2010-09-13

6 Comentarios

  1. 30

    Mostrar a su jefe el siguiente :

    set.seed(100)
    x1 <- runif(100,0,1)
    x2 <- as.factor(sample(letters[1:3],100,replace=T))
    
    y <- x1+x1*(x2=="a")+2*(x2=="b")+rnorm(100)
    summary(lm(y~x1*x2))

    Que da :

                Estimate Std. Error t value Pr(>|t|)    
    (Intercept)  -0.1525     0.3066  -0.498  0.61995    
    x1            1.8693     0.6045   3.092  0.00261 ** 
    x2b           2.5149     0.4334   5.802 8.77e-08 ***
    x2c           0.3089     0.4475   0.690  0.49180    
    x1:x2b       -1.1239     0.8022  -1.401  0.16451    
    x1:x2c       -1.0497     0.7873  -1.333  0.18566    
    ---
    Signif. codes:  0***0.001**0.01*0.05 ‘.’ 0.1 ‘ ’ 1 

    Ahora, con base en los valores de p se excluiría que uno? x2 es más significativa y no significativa en el mismo tiempo.


    Editar : Para aclarar : Este exaxmple no es la mejor, como se indica en los comentarios. El procedimiento en Stata y SPSS es AFAIK no esta basado en los valores de p de la prueba de T de los coeficientes, pero en el F-test después de la eliminación de una de las variables.

    Tengo una función que hace exactamente eso. Esta es una selección en el «p-valor», pero no de la prueba de T de los coeficientes o en el análisis de varianza de los resultados. Bien, siéntase libre de usar si se ve útil para usted.

    #####################################
    # Automated model selection
    # Author      : Joris Meys
    # version     : 0.2
    # date        : 12/01/09
    #####################################
    #CHANGE LOG
    # 0.2   : check for empty scopevar vector
    #####################################
    # Function has.interaction checks whether x is part of a term in terms
    # terms is a vector with names of terms from a model
    has.interaction <- function(x,terms){
    out <- sapply(terms,function(i){
    sum(1-(strsplit(x,":")[[1]] %in% strsplit(i,":")[[1]]))==0
    })
    return(sum(out)>0)
    }
    # Function Model.select
    # model is the lm object of the full model
    # keep is a list of model terms to keep in the model at all times
    # sig gives the significance for removal of a variable. Can be 0.1 too (see SPSS)
    # verbose=T gives the F-tests, dropped var and resulting model after 
    model.select <- function(model,keep,sig=0.05,verbose=F){
    counter=1
    # check input
    if(!is(model,"lm")) stop(paste(deparse(substitute(model)),"is not an lm object\n"))
    # calculate scope for drop1 function
    terms <- attr(model$terms,"term.labels")
    if(missing(keep)){ # set scopevars to all terms
    scopevars <- terms
    } else{            # select the scopevars if keep is used
    index <- match(keep,terms)
    # check if all is specified correctly
    if(sum(is.na(index))>0){
    novar <- keep[is.na(index)]
    warning(paste(
    c(novar,"cannot be found in the model",
    "\nThese terms are ignored in the model selection."),
    collapse=" "))
    index <- as.vector(na.omit(index))
    }
    scopevars <- terms[-index]
    }
    # Backward model selection : 
    while(T){
    # extract the test statistics from drop.
    test <- drop1(model, scope=scopevars,test="F")
    if(verbose){
    cat("-------------STEP ",counter,"-------------\n",
    "The drop statistics : \n")
    print(test)
    }
    pval <- test[,dim(test)[2]]
    names(pval) <- rownames(test)
    pval <- sort(pval,decreasing=T)
    if(sum(is.na(pval))>0) stop(paste("Model",
    deparse(substitute(model)),"is invalid. Check if all coefficients are estimated."))
    # check if all significant
    if(pval[1]<sig) break # stops the loop if all remaining vars are sign.
    # select var to drop
    i=1
    while(T){
    dropvar <- names(pval)[i]
    check.terms <- terms[-match(dropvar,terms)]
    x <- has.interaction(dropvar,check.terms)
    if(x){i=i+1;next} else {break}              
    } # end while(T) drop var
    if(pval[i]<sig) break # stops the loop if var to remove is significant
    if(verbose){
    cat("\n--------\nTerm dropped in step",counter,":",dropvar,"\n--------\n\n")              
    }
    #update terms, scopevars and model
    scopevars <- scopevars[-match(dropvar,scopevars)]
    terms <- terms[-match(dropvar,terms)]
    formul <- as.formula(paste(".~.-",dropvar))
    model <- update(model,formul)
    if(length(scopevars)==0) {
    warning("All variables are thrown out of the model.\n",
    "No model could be specified.")
    return()
    }
    counter=counter+1
    } # end while(T) main loop
    return(model)
    }
    • Buen ejemplo – podría haber utilizado este hace un par de meses!
    • Creo que es obvio que la decisión no está basada exclusivamente en los p-valores. Empezar por la eliminación de los menos significativos de orden más alto de las interacciones y, a continuación, el potencial cuadrático términos, antes de considerar cualquiera de las variables explicativas para la eliminación.
    • supuesto que no lo es. Pero quitar el término de interacción y verás que la situación sigue siendo la misma. Así que en realidad debería ir a un anova, pero luego te encuentras con el problema de que tipo. Y eso es una lata de gusanos que no voy a abrir.
    • No hay posibilidad argueing con el jefe en un nivel que suffisticated acerca de las estadísticas 😉 sin Embargo, voy a tratar de hacer la selección paso a paso «a mano» hasta entonces. Es raro igual que no hay una función de R para que, aun cuando el enfoque utilizando los valores de p viene desde el momento en que no eran demasiado altos costes computacionales de computiong la AIC …
    • Yo sé, está lejos de ser el mejor ejemplo. Era sólo para llevar de un punto a través. En cualquier caso, el «p-valor» método se basa en un test F, no en la prueba de t de los coeficientes. Por lo que la función de agregado debe dar el OP lo que él quiere.
    • Gracias por la gran ayuda. Tienes razón, para el «p-valor» método debo utilizar la Prueba F, se olvidó de eso. Gracias también por la publicación de su código, que también ayuda mucho.
    • De gran ayuda. Muchas gracias…
    • Meys hay una manera de limitar el número de variables ? por ejemplo, si desea seleccionar sólo 5 o 10 , ¿cómo puedo hacerlo?

  2. 18

    ¿Por qué no tratar de usar el step() función de especificar el método de prueba?

    Por ejemplo, para la eliminación hacia atrás, se escribe sólo un comando:

    step(FullModel, direction = "backward", test = "F")

    y para la selección paso a paso, simplemente:

    step(FullModel, direction = "both", test = "F")

    Este puede mostrar la AIC, así como valores de la F y los valores de P.

    • Tenga en cuenta que este procedimiento va a agregar/quitar la variable más importante de la prueba de valor. Pero la decisión de continuar o detener sigue basándose en la AIC. Creo que no es posible detener el proceso una vez que no hay variables significativas (como el OP quiere).
  3. 10

    Aquí es un ejemplo. A empezar con el más complicado de modelo: esto incluye la interacción entre las tres variables explicativas.

    model1 <-lm (ozone~temp*wind*rad)
    summary(model1)
    Coefficients:
    Estimate Std.Error t value Pr(>t)
    (Intercept) 5.683e+02 2.073e+02 2.741 0.00725 **
    temp          -1.076e+01 4.303e+00 -2.501 0.01401 *
    wind          -3.237e+01 1.173e+01 -2.760 0.00687 **
    rad           -3.117e-01 5.585e-01 -0.558 0.57799
    temp:wind      2.377e-01 1.367e-01 1.739 0.08519   
    temp:rad       8.402e-03 7.512e-03 1.119 0.26602
    wind:rad       2.054e-02 4.892e-02 0.420 0.47552
    temp:wind:rad -4.324e-04 6.595e-04 -0.656 0.51358

    La interacción de tres vías está claro que no es significativo. Esta es la forma de eliminarlo, para comenzar el proceso de modelo de simplificación:

    model2 <- update(model1,~. - temp:wind:rad)
    summary(model2)

    Dependiendo de los resultados, se puede continuar con la simplificación de su modelo:

    model3 <- update(model2,~. - temp:rad)
    summary(model3)
    ...

    Alternativamente, usted puede utilizar el modelo automático de simplificación de la función step, a ver
    qué bien hace:

    model_step <- step(model1)
    • Gracias por tu respuesta. Me acabo de dar cuenta que no he estado a mi pregunta clearlly suficiente: estoy buscando un comando o un paquete que hace esto automáticamente, con el p-valor como criterio. El comando de «paso» de los usos de la AIC. Yo también podría hacerlo «a mano» como usted suggets o programa, algunos de rutina que hace esto automáticamente, pero ya implementado la versión de que habría de venir en muy práctico.
  4. 7

    Si usted está tratando de conseguir el mejor modelo de predicción, entonces tal vez no importa demasiado, pero para nada más, no te molestes con este tipo de selección de modelo. Es equivocado.

    Utilizar un encogimiento de métodos tales como la regresión ridge (en lm.ridge() en el paquete de la MASA, por ejemplo), o el lazo, o la elasticnet (una combinación de cadena y lazo restricciones). De estos, sólo el lazo y red elástica a hacer algún tipo de selección de modelo, es decir, la fuerza de los coeficientes de algunas covariables a cero.

    Ver la Regularización y la Contracción de la sección de la La Máquina De Aprendizaje vista de tareas en CRAN.

  5. 0

    Como se ha mencionado por Gavin Simpson la función fastbw de rms paquete puede ser usado para seleccionar las variables con el p-valor. Abajo es un ejemplo utilizando el ejemplo dado por George Dontas. Utilice la opción rule='p' para seleccionar el valor p de los criterios.

    require(rms)
    model1 <- ols(Ozone ~ Temp * Wind * Solar.R, data=airquality)
    model2 <- fastbw(fit=model1, rule="p", sls=0.05)
    model2
    • En el texto, que significó la rms el paquete. También hay un rsm el paquete, por lo que esto podría causar confusión.
    • Gracias rvl por la aclaración. Me fijo.

Dejar respuesta

Please enter your comment!
Please enter your name here