[R] AIC function and Step function

Dana77 luckyinwind at yahoo.com
Mon Dec 1 01:05:35 CET 2008


 Thanks for kind help from Steven and Christos last time.  Now I got new
problem regarding the codes for calculating the "weights" (w) in "AIC ()
function". 
The original code is as below:
 > getAnywhere("logLik.lm")
function (object, REML = FALSE, ...) 
  {
    res <- object$residuals
    p <- object$rank
    N <- length(res)
    if (is.null(w <- object$weights)) {
        w <- rep.int(1, N)
    }    else {
        excl <- w == 0
        if (any(excl)) {
            res <- res[!excl]
            N <- length(res)
            w <- w[!excl]
        }
    }

 Now my question is, if I use "lm()" function to fit a multiple linear
regression model, such as "mod.fit<-lm(formula = Y~ X1 + X2 + X3, data =
set1)", what code could I use to extract the "weights" (w) out? or how to
calculate the weights(w) shown in above codes? Thanks for your time and kind
help!

Dana



Steven McKinney wrote:
> 
> Hi Dana,
> 
> Many thanks to Christos Hatzis who sent
> me an offline response, pointing out the
> new functions that make this much
> easier than my last suggestions:
> methods() and getAnywhere()
> 
>> methods("extractAIC")
> [1] extractAIC.aov*     extractAIC.coxph*   extractAIC.glm*    
> extractAIC.lm*      extractAIC.negbin* 
> [6] extractAIC.survreg*
> 
>    Non-visible functions are asterisked
>> getAnywhere("extractAIC.coxph")
> A single object matching ‘extractAIC.coxph’ was found
> It was found in the following places
>   registered S3 method for extractAIC from namespace stats
>   namespace:stats
> with value
> 
> function (fit, scale, k = 2, ...) 
> {
>     edf <- length(fit$coef)
>     loglik <- fit$loglik[length(fit$loglik)]
>     c(edf, -2 * loglik + k * edf)
> }
> <environment: namespace:stats>
>> 
> 
> Thank you Christos.
> 
> 
> That said, one of the advantages of getting
> the source code is that it has comments that
> are stripped out when the code is sourced into R
> 
> e.g. from the command line
> 
>> getAnywhere(AIC.default)
> A single object matching ‘AIC.default’ was found
> It was found in the following places
>   registered S3 method for AIC from namespace stats
>   namespace:stats
> with value
> 
> function (object, ..., k = 2) 
> {
>     ll <- if ("stats4" %in% loadedNamespaces()) 
>         stats4:::logLik
>     else logLik
>     if (length(list(...))) {
>         object <- list(object, ...)
>         val <- lapply(object, ll)
>         val <- as.data.frame(t(sapply(val, function(el) c(attr(el, 
>             "df"), AIC(el, k = k)))))
>         names(val) <- c("df", "AIC")
>         Call <- match.call()
>         Call$k <- NULL
>         row.names(val) <- as.character(Call[-1])
>         val
>     }
>     else AIC(ll(object), k = k)
> }
> <environment: namespace:stats>
> 
>>From the source file
> 
> 
> #  File src/library/stats/R/AIC.R
> #  Part of the R package, http://www.R-project.org
> #
> #  This program is free software; you can redistribute it and/or modify
> #  it under the terms of the GNU General Public License as published by
> #  the Free Software Foundation; either version 2 of the License, or
> #  (at your option) any later version.
> #
> #  This program is distributed in the hope that it will be useful,
> #  but WITHOUT ANY WARRANTY; without even the implied warranty of
> #  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
> #  GNU General Public License for more details.
> #
> #  A copy of the GNU General Public License is available at
> #  http://www.r-project.org/Licenses/
> 
> #### Return the object's value of the Akaike Information Criterion
> #### (or "An Inf.. Crit..")
> 
> AIC <- function(object, ..., k = 2) UseMethod("AIC")
> 
> ## AIC for logLik objects
> AIC.logLik <- function(object, ..., k = 2)
>     -2 * c(object) + k * attr(object, "df")
> 
> AIC.default <- function(object, ..., k = 2)
> {
>     ## AIC for various fitted objects --- any for which there's a logLik()
> method:
>     ll <- if("stats4" %in% loadedNamespaces()) stats4:::logLik else logLik
>     if(length(list(...))) {# several objects: produce data.frame
> 	object <- list(object, ...)
> 	val <- lapply(object, ll)
> 	val <- as.data.frame(t(sapply(val,
> 				      function(el)
> 				      c(attr(el, "df"), AIC(el, k = k)))))
> 	names(val) <- c("df", "AIC")
>         Call <- match.call()
>         Call$k <- NULL
> 	row.names(val) <- as.character(Call[-1])
> 	val
>     } else AIC(ll(object), k = k)
> }
> 
> 
> 
> Steven McKinney
> 
> Statistician
> Molecular Oncology and Breast Cancer Program
> British Columbia Cancer Research Centre
> 
> email: smckinney +at+ bccrc +dot+ ca
> 
> tel: 604-675-8000 x7561
> 
> BCCRC
> Molecular Oncology
> 675 West 10th Ave, Floor 4
> Vancouver B.C. 
> V5Z 1L3
> Canada
> 
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide
> http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
> 
> 

-- 
View this message in context: http://www.nabble.com/AIC-function-and-Step-function-tp20728043p20764062.html
Sent from the R help mailing list archive at Nabble.com.



More information about the R-help mailing list