L.S.
After moment fitting, which I programmed as follows:
#Function which, given the first two moment, computates the density,
#according to momentfitting
momentfit<-function(ES,c){
res<-(1:5)
if(c<=1){
k<-floor(1/c^2)
p<-((1)/(1+c^2))*((k+1)*c^(2)-sqrt((k+1)*(1+c^(2))-(k+1)^(2)*c^(2)))
mu<-(k+1-p)/ES
res[1]<-0
res[2]<-k
res[3]<-p
res[4]<-mu
res[5]<-0
}
if(c>1){
p1<-((1)/(2))*(1+sqrt((c^2-1)/(c^2+1)))
p2<-1-p1
mu1<-2*p1/ES
mu2<-2*p2/ES
res[1]<-1
res[2]<-p1
res[3]<-p2
res[4]<-mu1
res[5]<-mu2
}
return(res)
}
I have the either one of the following distributionfunctions:
ES<-a
c<-b
para<-momenfit(a,c)
if(para[1]==0){
fun1<-function(x){(para[3]*dgamma(x,para[2],para[4])+(1-para[3])*dgamma(x,para[2]+1,para[4]))}
}
if(para[1]==1){
fun2<-function(x){(para[2]*dgamma(x,1,para[4])+para[3]*dgamma(x,1,para[5]))}
}
I now want to generate a random number according to the above distrubutions. Can
somebody help me? Thanks.
Martin