An embedded and charset-unspecified text was scrubbed... Name: not available URL: <https://stat.ethz.ch/pipermail/r-help/attachments/20111112/72f1912c/attachment.pl>
Please Help
2 messages · avish kamdar, David Winsemius
Two problems with your post: a) you used a format that is not accepted by the mail server and resulted in all of your linefeeds being removed, so it is now basically unreadable. As the Posting Guide says, please use plain text. (And please read the Posting Guide for yourself ... and do it now before reading further.) All the free mailing services provide mechanisms for constructing plain text emails. b) it appears to be homework, and homework submission are discouraged. You can attach a plea describing why you should be exempted from this rule, but it might require identifying your self more completely and explaining why you are not using the services of the educational institution you are studying under.
On Nov 12, 2011, at 8:09 AM, avish kamdar wrote:
HiI want to construct a logliikelood function in RHere is the
situationy=number of particles emitted in 1 hr
period~pois(30)p=probability of detection of radiation
particlesx=number of particles detected by a radiation
detector~pois(30p)where p~beta(a,1)I have to calculate the
loglikehood for a for the range a(2,50)I wish to simulate 100 random
samples for each aHere is my code:-m=481n=100x = c(15, 36, 29, 28,
37, 32, 25, 27, 31, 21, 25, 27, 28, 31, 28, 20, 34, 25, 20, 34,
15,21, 28, 24, 31, 19, 34, 29, 18, 25, 16, 19, 44, 26, 34, 31, 21,
28, 11, 31, 21, 34, 25, 25,30, 23, 21, 35, 36, 21, 27, 29, 30, 22,
25, 30, 24, 27, 28, 22, 36, 29, 33, 35, 30, 32, 27,26, 25, 27, 23,
21, 39, 33, 24, 21, 19, 34, 32, 28, 27, 28, 23, 20, 24, 29, 21, 22,
31, 28,27, 28, 29, 21, 30, 28, 31, 22, 29, 18)a=c(1:m)*0.1+1.9p<-
rbeta(n,a,1)loglik<-numeric(m)for(i in 1:n){p[i]<-rbeta(n,a[i],
1)x[i]<-rpois(n,p[i]*30)loglik<-function(x,p)logLik(a=a(i))<-n*log(a)
+(log(30)*sum(x))-sum(log(factorial(x)))-30*(!
sum(p[i]))+sum(x[i]*log(p[i]))+(a[j]-1)*sum(log(p[i]))}But somehow I
am not getting the answer please HELPP.S. the data for x is a random
sample of number of particles detected in 1 hr (given in the
question)Really appreciate the helpThanks
[[alternative HTML version deleted]]
David Winsemius, MD West Hartford, CT