Jacob Warren (RIT Student)
2014-Apr-24 19:27 UTC
[R] LogLikelihood of a Distribution Given Fixed Parameters
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(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) [[alternative HTML version deleted]]
Rolf Turner
2014-Apr-25 04:50 UTC
[R] LogLikelihood of a Distribution Given Fixed Parameters
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(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)