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.
if (any(excl)) {
res <- res[!excl]
N <- length(res)
w <- w[!excl] }
}
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))))
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.
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 readthe> code, I found there are two cases for "w(weights)". If "w" is notzero,> 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 linesafter 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 notzero> 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
Thank you, Steven. It helps! Best, Dana Steven McKinney wrote:> > 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 > > ______________________________________________ > 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/Help-with-reading-code-tp20823979p20847549.html Sent from the R help mailing list archive at Nabble.com.