AIC function and Step function

On Sun, Nov 30, 2008 at 5:05 PM, Dana77 <luckyinwind at yahoo.com> wrote:

 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?

mod.fit won't have weights because you didn't specify any through the
weights argument to lm.  If you had, you could extract them using the
same technique used in the above code:  w <- mod.fit$weights

hth,

Kingsford Jones

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

--
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.