search for: lgamma

Dear R users, i try to integrate lgamma from 0 to Inf. But here i get the message "roundoff error is detected in the extrapolation table", if i use 1.0e120 instead of Inf the computation works, but this is against the suggestion of integrates help information to use Inf explicitly. Using stirlings approximation doesnt bring the...

...975 then ucl(m)= y-1 how do I code these instructions into this code. thanks Aruike pp=function(x,n,M){z=1.0;a=2.3071430;b=7.266064;H=3 out.h=c() out.y=c() out.m=c() out.prob=c() for(h in 1:H){ for(m in 1:M){ for(y in 0:m){ g=lgamma(m+z)+lgamma(n[h]+a+b)+lgamma(x[h]+y+a)+lgamma(n[h]+m+b-x[h]-y) g=g-lgamma(y+z)-lgamma(m-y+z)-lgamma(x[h]+a)-lgamma(n[h]+b-x[h])- lgamma(n[h]+m+a+b) out.h=c(out.h,h) out.y=c(out.y,y) out.m=c(out.m,m) out.prob=c(out.prob, exp(g)) out.cum=c(cumsum(out.prob))...

...compile but I haven't given up hope on that one yet) -- the last compile I did was 0.61.1, which worked smoothly. Running ./configure --prefix=/usr/local/apps/R/R-0.62.3 [--with-g77] (tried both with and without g77, eventually get to the same place), I get a conflict with the function lgamma being defined both in R and in the system's libm.a file: make [snip] cd main; make g77 -o R.binary arithmetic.o [snip] version.o ../lib/libgraphics.a ../lib/libunix.a ../lib/libappl.a ../lib/libmath.a -lSM -lICE -L/usr/lib/X11 -lX11 -ldl -ltermcap -lm collect2: ld returned 2 exit status...

I am porting some FORTRAN to R in which an Inf triggers an if(). The trigger is infinite on exp(lgamma(OVER)). What is the canonical R style of determining OVER when exp(OVER)== Inf? The code structure that I am porting is best left intact--so I need to query R somehow to the value of OVER that causes exp(lgamma(OVER)) to equal Inf. On my system, exp(lgamma(171)) is about first to equal Inf. I...

....000 but the 2nd and 3rd m returned 2.000 and 3.000 instead of 1.000. thanks Aruike pp=function(x,n,M){z=1.0;a=2.3071430;b=7.266064;H=3 out.h=c() out.y=c() out.m=c() out.prob=c() for(h in 1:H){ for(m in 1:M){ for(y in 0:m){ g=lgamma(m+z)+lgamma(n[h]+a+b)+lgamma(x[h]+y+a)+lgamma(n[h]+m+b-x[h]-y) g=g-lgamma(y+z)-lgamma(m-y+z)-lgamma(x[h]+a)-lgamma(n[h]+b-x[h])- lgamma(n[h]+m+a+b) out.h=c(out.h,h) out.y=c(out.y,y) out.m=c(out.m,m) out.prob=c(out.prob, exp(g)) out.cum=c(cumsum(out.prob))...

...: $ operator is invalid for atomic vectors In addition: There were 11 warnings (use warnings() to see them) > warnings() Warning messages: 1: In optim(param, fn = mt.dev, method = method, control = control, ... : unknown names in control: tol 2: In fn(par, ...) : value out of range in 'lgamma' 3: In fn(par, ...) : value out of range in 'lgamma' 4: In fn(par, ...) : value out of range in 'lgamma' 5: In fn(par, ...) : value out of range in 'lgamma' 6: In fn(par, ...) : value out of range in 'lgamma' 7: In fn(par, ...) : value out of range in 'lgamma...

.... It is very appreciate if you can help. function (parameters,y,x1,x2) { p<-parameters[1] alpha1<-parameters[2] beta1<-parameters[3)] delta1<-parameters[4] alpha2<-parameters[5] mu<-alpha1*((x1)^beta1)*exp(-delta1*(x1^alpha2)) if(y>0 & x1>0 & x2==1, L<-lgamma(y+p)+p*(log(p)-log(mu+p))+y*(log(mu)-log(mu+p))-lfactorial(y)-lgamma(p) ) else if(y>0 & x1>0 & x2==2, L<-lgamma(y+p)+p*(log(p)-log(mu+p))+y*(log(mu)-log(mu+p))-lfactorial(y)-lgamma(p) ) else if(y>0 & x1>0 & x2==3, L<-lgamma(y+p)+p*(log(p)-log(mu+p))+y*(log(...

Hi all, why factorial(150) shows the error out of range in 'gammafn'? I have to calculate the number of subset formed by 150 samples taking 10 at a time. How is this possible? best [[alternative HTML version deleted]]

Hello, I submitted a package which used Fortran functions isnan and lgamma. However, I was told that: isnan and lgamma are not Fortran 95 functions. I was asked to write 'cross-platform portable code' and so should not be writing GNU extensions to Fortran. See http://cran.r-project.org/web/checks/check_results_mpath.html, which will shortly show installation fa...

The following seems very strange: > lgamma(-0.04) [1] 3.243307 > gamma(-0.04) [1] -25.61830 Jim -.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.- r-devel mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html Send "info", "help", or "[un]subscribe" (in the &quot...

...milar expression as te left hand except k becomes k+1 my code is below, I don't know why R keep telling me the syntax wrong but I am not clear, thanks. f<-function(u,k){(1+u^2/(k-1))^(-k/2)} g<-function(u,k){(1+u^2/k)^(-(k+1)/2)} for(k in 4:25){ out<-uniroot(function(a){ m<-exp(lgamma(k/2)-lgamma((k-1)/2))*1/sqrt(k-1) n<-exp(lgamma((k+1)/2)-lgamma(k/2)*1/sqrt(k) m*integrate(f,lower=0,upper=sqrt(a^2*(k-1)/(k-a^2)),rel.tol=.Machine$double.eps^0.25) -n*integrate(g,lower=0,upper=sqrt(a^2*k/(k+1-a^2)),rel.tol=.Machine$double.eps^0.25) } lower=1.4, upper=2) print(unlist(out))} --...

Dear R users: I wonder anyone is aware of such a thing from lgamma() > lgamma(-0.8) [1] 1.747207 I thought it should be NA as in S-PLUS. Both R-1.3.1 and R-1.5.1 report this result. I don't have the latest R-1.6.1 so I don't know whether this is corrected or not. Happy holidays, Paul.

...e instructions into this code. thanks Aruike pp=function(x,n,M){z=1.0;a=2.3071430;b=7.266064;H=3 out.h=c() out.y=c() out.m=c() out.prob=c() for(h in 1:H){ for(m in 1:M){ for(y in 0:m){ g=lgamma(m+z)+lgamma(n[h]+a+b)+lgamma(x[h]+y+a)+lgamma(n[h]+m+b-x[h]-y) g=g-lgamma(y+z)-lgamma(m-y+z)-lgamma(x[h]+a)-lgamma(n[h]+b-x[h])- lgamma(n[h]+m+a+b) out.h=c(out.h,h) out.y=c(out.y,y) out.m=c(out.m,m) out.prob=c(out.prob, exp(g)) out.cum=c(cumsum(out....

Hi, I found the following on Windows 2000/NT R Version 1.9.1 (2004-06-21) (also Version 1.9.0): The S4 group "Math" doesn't work as documented; i.e., "log", "log10", "gamma" and "lgamma" are included in the documentation but don't work. See example code below. Moreover, what about 'genericForPrimitive' which is used in 'getGeneric'. It seems that this method is not included in the R Version 1.9.1 (also 1.9.0). See the example code of John Chambers at the...

...program. It works for me. Tthe program is described as below: function (parameters,y1,x11) { p<-parameters[1] alpha1<-parameters[2] beta1<-parameters[3] delta1<-parameters[4] lamda1<-parameters[5] mu<-alpha1*((x11)^beta1)*exp(-delta1*(x11^lamda1)) ifelse(y1>0|x11>0, L<-lgamma(y1+p)+p*(log(p)-log(mu+p))+y1*(log(mu)-log(mu+p))-lfactorial(y1)-lgamma(p) ,Inf) L } This is working for me. 2 Way two: I select the data whose x2=1 in the whole range of data. It works but it is not right comparing the value of MLE. the program is: function (parameters,y,x1,x2) { p&l...

...the method described in the Guo paper is available at http://www.stat.unipg.it/stat/statlib/general/negmul. I''ve been approaching this problem by modifying the loglik glm.nb function as such: ## original loglik function in glm.nb ## loglik <- function(n, th, mu, y) { ## sum(lgamma(th + y) - lgamma(th) - lgamma(y + 1) + th * ## log(th) + y * log(mu + (y == 0)) - (th + y) * log(th + ## mu)) ## } loglik <- function(n, th, mu, y, mu.grouped, y.grouped) { sum( th*log(th))*length(y.grouped) + sum(y*log(mu+ (y == 0)))+ sum(lgamma(y.groupe...

I have to calculate some formula like: gamma(x)/(gamma(x+y) and I observed that for relatively big values of x, R returns infinity and so cannot compute the formula. Is it possible to force R to give the real value of gamma(x) instead of Inf ? thanks

Dear R-users, Is there an R function to compute the multinomial beta function? That is, the normalizing constant that arises in a Dirichlet distribution. For example, with three parameters the beta function is Beta(n1,n2,n2) = Gamma(n1)*Gamma(n2)*Gamma(n3)/Gamma(n1+n2+n3) Thanks in advance for any assisstance. Regards, Greg [[alternative HTML version deleted]]

Dear R users, I have to calculate gamma functions for negative numbers beyond -171.4. e.x. gamma(-500.4) I got following: > gamma(-170.4) [1] -5.824625e-308 > gamma(-171.4) [1] 0 Warning message: underflow occurred in 'gammafn' I have tried to use a recursion getting values a little futher -180. How could I solve this problem? Thank you beforehand. Chuse.

...z tell me how to avoid this "for" loop blow?? ---------------------------------------------------------------------------------------------- for (k in 1:n){ ymax <- max( y1[k], y2[k] ) i <- 0:ymax sums<- -lgamma(y1[k]-i+1)-lgamma(i+1)-lgamma(y2[k]-i+1) maxsums <- max(sums) sums <- sums - maxsums lsum <- log( sum(exp(sums)) ) + maxsums lbp[k] <- y1[k]*log(4) + y2[k]*log(5) + lsum } elbp <- exp(...