4
$\begingroup$

I'm trying to fit the GAMLSS library's Sichel distribution to some zero-truncated data, but the only way to get the function to work is to include the zero-class anyway but give it a frequency of 0, which doesn't take into account the zero-truncated nature of my data. Can anyone suggest a way to properly "redistribute" the zero-class's probability to the remaining probabilities (or some other, better, course of action using Sichel)?

If you run the following example, you'll see that sum(pdf2) equals 1, but that the zero class that has a probability in my case of 0 is still allocated around 27% of the cum probability: 

Counts = data.frame(n = c(0,1,2,3,4,5,6,7,8,9,10),
                    freq = c(0,182479,76986,44859,24315,49,100,490,106,0,2))

gamlss(n~1,family=SICHEL, control=gamlss.control(n.cyc=50),data=Counts )

pdf2 = dSICHEL(x=with(Counts, n), mu = 1.610, sigma = 98.43, nu = 3.315)

print( with(Counts, cbind(n, freq, fitted=pdf2*sum(freq))), dig=9)

sum(pdf2)
Improve this question
$\endgroup$

3 Answers 3

Reset to default
2
$\begingroup$

In R, load the package gamlss.tr (for fitting truncated distributions), then: 

y<-rSICHELtr(1000,2,0.5,-1)
hist(y)
gen.trun(par=0,family="SICHEL", type="left",delta=0.0001)
m1<-gamlss(y~1,family=SICHELtr)
summary(m1)

The above code generates a sample of 1000 from the truncated SICHEL distribution, then fits it. Note in the simulation mu=2, sigma=0.5 and nu=-1.

In the summary fitted model gives fitted log(mu), log(sigma) and nu.

Improve this answer
$\endgroup$
1
  • $\begingroup$ Welcome to the site, @Robert. We prefer if you don't sign your posts; note that your avatar & username (w/ a link to your userpage) are added automatically to your posts. Since you're new here, you'll want to read our about page & our FAQ, which contain info like this about CV. $\endgroup$
    – gung - Reinstate Monica
    Commented Mar 3, 2013 at 18:07
1
$\begingroup$

By zero-truncated, do you mean that any data that would have had a 0 as a response is just missing? In that case, can't you just put the 0s in?

Or do you mean that some proportion of the time, instead of getting a sensical answer, you get a 0 instead? That sounds like zero-inflation to me. In that case, there are zero-inflated poisson and similar in GAMLSS.

I don't know of a zero-inflated Sichel, and there's nothing in GAMLSS for it, but is there a particularly good reason for using the Sichel distribution for your data? Does it reflect the underlying process particularly well? (I believe that the Sichel represents a mixture model of Poissons, where the meta-distribution is distributed as Inverse Gaussian...)

Improve this answer
$\endgroup$
4
  • $\begingroup$ Follow-up: 1. The SICHEL() function is in the gamlss.dist package. 2. Would a simple (and correct!) method to adjust the > 0 probabilities be to divide them by (1-.27)? --Jeff $\endgroup$
    – Jeff Tyzzer
    Commented Jan 11, 2011 at 0:11
  • $\begingroup$ By zero-truncated, I mean that the probability of a zero-count is zero. My data are the count of items purchased per transaction at a convenience store over a five-day period--obviously there's at least one item being bought per purchase (but usually fewer than, say, 20). The schema looks something like: DATE, TRANSACTION_ID, COUNT_OF_PURCHASES As to Sichel, I do have reason to believe that that distribution might provide a good fit, but I'm trying several others as well (e.g., ZTNB). Thanks for your reply. --Jeff $\endgroup$
    – Jeff Tyzzer
    Commented Jan 26, 2011 at 20:38
  • $\begingroup$ Oh, I see. ZTNB seems viable. Or you could model visits including unobserved non-purchases, as you originally noted. That would basically be a free parameter, but if Sichel is describing some other aspect of the underlying behavior, it might be worth it. $\endgroup$
    – Harlan
    Commented Jan 26, 2011 at 20:45
  • $\begingroup$ Incidentally, you should have responded to my answer, not in a separate answer. And if you like the answer, please up-vote it, and if it actually answers the question to your satisfaction, accept it. (Likewise, your follow-up should have been a comment on your question, not an answer.) Welcome to Stack Exchange! :) $\endgroup$
    – Harlan
    Commented Jan 26, 2011 at 20:46
1
$\begingroup$

The distribution of the zero truncated law is $$ \mathbb P(X = k | X > 0) = { \mathbb P(X = k) \over 1 - \mathbb P(X=0) },$$ for $k>0$.

So in R it is given by dSICHEL( x, mu, sigma, nu)/(1-dSICHEL(0, mu, sigma, nu)) and its log by dSICHEL( x, mu, sigma, nu, log=TRUE)-log(1-dSICHEL(0, mu, sigma,nu))... it is then easy to compute the log-likelihood for given parameters $\mu, \sigma, \nu$. You could try to maximize it using nlm or optim.

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.