[R] Help with reading code

[Steven McKinney]

Hi Dana

> -----Original Message-----
> From: r-help-bounces at r-project.org
[mailto:r-help-bounces at r-project.org]
> On Behalf Of Dana77
> Sent: Wednesday, December 03, 2008 3:24 PM
> To: r-help at r-project.org
> Subject: [R] Help with reading code
> 
> 
> I would like to give out the equation for calculating the maximum
> likelihood.
> Below is the code, but I still have problems with it.  After I read
the
> code, I found there are two cases for "w(weights)".  If  "w" is not
zero,
> then the equation is given as "val <- 0.5 * (sum(log(w)) - N * (log(2
*
> pi)
> + 1 - log(N) +
>         log(sum(w * res^2))))". However, if "w" is zero, then I do not
> know
> what equation it should be since it does not make any sense for
"log0".
> Hope
> someone can help me to figure this out. Thanks!
> 
> 
> 
> 
> 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  #####I can not understand the following lines
after this.

The command above sets a variable to "exclude" with.
Any observation with weight equal to zero will get excluded.
excl will have value TRUE for all observations with weight 0.

>         if (any(excl)) { ### If there are any observations to exclude
>             res <- res[!excl]   ### then keep the ones with weight not
zero
>             N <- length(res)
>             w <- w[!excl]        }

so the above chunk keeps only the residuals and weights
for the observations with non-zero weights

>     }
>     N0 <- N
>     if (REML)
>         N <- N - p
>     val <- 0.5 * (sum(log(w)) - N * (log(2 * pi) + 1 - log(N) +
>         log(sum(w * res^2))))

so now there are no more observations with w == 0 in the above equation

>     if (REML)
>         val <- val - sum(log(abs(diag(object$qr$qr)[1:p])))
>     attr(val, "nall") <- N0
>     attr(val, "nobs") <- N
>     attr(val, "df") <- p + 1
>     class(val) <- "logLik"
>     val
> }
> 
> 
> Best,
> 
> Dana
> --
> View this message in context: http://www.nabble.com/Help-with-reading-
> code-tp20823979p20823979.html
> Sent from the R help mailing list archive at Nabble.com.
> 
> ______________________________________________
> 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.



 

Steven McKinney

Statistician
Molecular Oncology and Breast Cancer Program
British Columbia Cancer Research Centre

email: smckinney at bccrc.ca
tel: 604-675-8000 x7561

BCCRC
Molecular Oncology
675 West 10th Ave, Floor 4
Vancouver B.C. 
V5Z 1L3

Canada