Indeed it should be lower=FALSE. Duhhhh! Sorry 'bout that! cheers, Rolf On 29/04/14 02:51, Jacob Warren (RIT Student) wrote:
Thanks Rolf. That took care of it. It should be lower=FALSE through right? I want the upper tail because my values are right censored? Regards, Jake On Fri, Apr 25, 2014 at 12:50 AM, Rolf Turner <[email protected] <mailto:[email protected]>> wrote: As usual I am too lazy to fight my way through the rather convoluted code presented, but it seems to me that you just want to calculate a log likelihood. And that is bog-simple: The log likelihood for i.i.d. data is just the sum of log f(y_i) where the y_i are your observed values and f() is the density function of the distribution that you have in mind. Where there is (right) censoring you take the sum of log f(y_i) over all the non-censored values and then add k*(1-F(cens.time)) where k is the number of censored values and F() is the cumulative distribution function corresponding to f(). In your case it would appear that f(y) = dlnorm(y,1.66,0.25) and F(y) = plnorm(y,1.66,0.25). Note that instead of using 1-F(cens.time) you can use plnorm(cens.time,1.66,0.25,__lower=TRUE) and that instead of taking logs explicitly you can set log=TRUE in the calls to dlnorm() and plnorm(). cheers, Rolf Turner On 25/04/14 07:27, Jacob Warren (RIT Student) wrote: I'm trying to figure out if there is a way in R to get the loglikelihood of a distribution fit to a set of data where the parameter values are fixed. For example, I want to simulate data from a given alternate lognormal distribution and then I will fit it to a lognormal distribution with null parameter values to see what the likelihood of the null distribution is given random data from the alternate distribution. I have been using fitdistrplus for other purposes but I cannot use it to fix both parameter values. Here is an example of what I've been working with... nullmu<-1.66 #set null mu altmu<-1.58 #set alt mu sd.log<-0.25 #set common sigma cens.time<-6 #if simulated times are greater than this turn them into right censored times #simulating lognormal data (time) from altnative dist (sim<-rlnorm(n=samplesize, meanlog=altmu, sdlog=sd.log)) #if the time was > cens.time replace time with cens.time (sim[which(sim>cens.time)]<-__cens.time) sim #create a variable indicating censoring (cens<-sim) cens[which(sim==cens.time)]<-__NA cens #create the data frame to be passed to fitdistcens and fitdist (x<-data.frame(left=sim,right=__cens)) #if there is censored data use fitdistcens else use fitdist ifelse(length(which(is.na <http://is.na>(__cens)))>0, simfit<-fitdistcens(censdata=__x, distr="lnorm"), simfit<-fitdist(data=x[,1], distr="lnorm") ) #Now I can get the loglikelihood of the MLE fitted distribution simfit$loglik #I want to get the loglikelihood of the distribution with the null parameterization #This is what I can't get to work #I can't seem to find any function that allows me to set both parameter values #so I can get to loglikelihood of the of the parameterization given the data nulldist<-fitdistcens(__censdata=x, distr="lnorm", start=list(meanlog=nullmu, sdlog=sd.log) #Then I want to do a likelihood ratio test between the two distributions pchisq((-2*simfit$loglik--2*__nulldist$loglik), df=2, lower.tail=FALSE)
______________________________________________ [email protected] 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.
