Re: [R] fitting truncated normal distribution

From: Sundar Dorai-Raj <sundar.dorai-raj_at_pdf.com>
Date: Fri 18 Aug 2006 - 00:47:18 EST

Hi, Markus,

Are these always integers? Why do you think they should be normal or Weibull? Seems more like a mixture with a point mass at 0 and something else (e.g. Poisson, negative binomial, normal). Though it's hard to tell with what you have provided. If that's the case you'll have to write your own likelihood function or, if they are integers, use zip (zero-inflated Poisson) or zinb (zero-inflated negative binomial). Do an RSiteSearch to find many packages will do these fits.

RSiteSearch("zero-inflated")

Again, this is pure speculation based on your "x" below alone and no other information (I'm not sure what "demand-data" means).

HTH, --sundar

Schweitzer, Markus wrote:
> Sorry, that I forgot an example.
>
> I have demand-data which is either 0 or a positive value.
>
> When I have an article which is not ordered very often, it could look
> like this:
>
> x=c(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1280,0,0,0,0,640,0,0
> ,0,0,0,0,0,0,0)
>
>
>

>>library(MASS) ## for fitdistr
>>library(msm) ## for dtnorm
>>
>>dtnorm0 <- function(x, mean = 0, sd = 1, log = FALSE) {
>>   dtnorm(x, mean, sd, 0, Inf, log)
>>}
>>fitdistr(x,dtnorm0,start=list(mean=0,sd=1))

>
>
> Unfortunately I get the same error message.
> I found a function, that works for a weibull distribution and tried to
> apply it but it didn't work neither
>
> # truncated weibull distribution
>
> #dweibull.trunc <-
> #function(x, shape, scale=1, trunc.=Inf, log=FALSE){
> # ln.dens <- (dweibull(x, shape, scale, log=TRUE)
> # -pweibull(trunc., shape, scale = 1, lower.tail = TRUE, log.p =
> #TRUE))
> # if(any(oops <- (x>trunc.)))
> # ln.dens[oops] <- (-Inf)
> # if(log)ln.dens else exp(ln.dens)
> #}
> #
> #x <- rweibull(100, 1)
> #range(x)
> #x4 <- x[x<=4]
> #fitdistr(x4, dweibull.trunc, start=list(shape=1, scale=1), trunc=4)
>
> ########################################################################
> ########
> # truncated normal distribution
>
> dtnorm0 <- function(x, mean, sd, a=0, log = FALSE) {
> ln.dens <- (dnorm(x, mean, sd)
> - pnorm(a, mean, sd, lower.tail=TRUE, log.p =TRUE))
>
> if(any(oops <- (x<a)))
> ln.dens[oops] <- (-Inf)
> if(log)ln.dens else exp(ln.dens)
> }
>
> fitdistr(x, dtnorm0, start = list(mean = 0, sd = 1))
>
> Maybe, when I alter mean and sd, I get an answer, which is not really
> satisfactory. I hope, there is a solution possible
> And thank you in advance
>
> markus
>
>
>
>
>
>
>
> Sorry, didn't notice that you *did* mention dtnorm is part of msm.
> Ignore that part of the advice...
>
> --sundar
>
> Sundar Dorai-Raj wrote:
>
>>aon.912182281.tmp@aon.at wrote:
>>
>>
>>>Hello,
>>>I am a new user of R and found the function dtnorm() in the package

>
> msm.
>
>>>My problem now is, that it is not possible for me to get the mean and

>
> sd out of a sample when I want a left-truncated normal distribution
> starting at "0".
>
>>>fitdistr(x,dtnorm, start=list(mean=0, sd=1))
>>>
>>>returns the error message 
>>>"Fehler in "[<-"(`*tmp*`, x >= lower & x <= upper, value = numeric(0))

>
> : nichts zu ersetzen"
>
>>>I don't know, where to enter the lower/upper value. Is there a

>
> possibility to program the dtnorm function by myself?
>
>>>Thank you very much in advance for your help, markus
>>>
>>>-------------------------------------------
>>>Versendet durch aonWebmail (webmail.aon.at)
>>>
>>>
>>>______________________________________________
>>>R-help@stat.math.ethz.ch 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.
>>
>>
>>Hi, Markus,
>>
>>You should always supply the package name where dtnorm is located. My 
>>guess is most don't know (as I didn't) it is part of the msm package.
>>Also, you should supply a reproducible example so others may 
>>understand your particular problem. For example, when I ran your code 
>>on data generated from "rtnorm" (also part of msm) I got warnings 
>>related to the NaNs generated in pnorm and qnorm, but no error as you 
>>reported. Both of these suggestions are in the posting guide (see

>
> signature above).
>
>>So, to answer your problem, here's a quick example.
>>
>>library(MASS) ## for fitdistr
>>library(msm) ## for dtnorm
>>
>>dtnorm0 <- function(x, mean = 0, sd = 1, log = FALSE) {
>>   dtnorm(x, mean, sd, 0, Inf, log)
>>}
>>
>>set.seed(1) ## to others may reproduce my results exactly x <- 
>>rtnorm(100, lower = 0) fitdistr(x, dtnorm0, start = list(mean = 0, sd 
>>= 1))
>>
>>Note, the help page ?fitdistr suggests additional parameters may be 
>>passed to the density function (i.e. dtnorm) or optim. However, this 
>>won't work here because "lower" is an argument for both functions. 
>>This is the reason for writing dtnorm0 which has neither a lower or an

>
>
>>upper argument.
>>
>>HTH,
>>
>>--sundar
>>
>>______________________________________________
>>R-help@stat.math.ethz.ch 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.

>
>
> ______________________________________________
> R-help@stat.math.ethz.ch 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.
>
>
>
>
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help@stat.math.ethz.ch 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.


R-help@stat.math.ethz.ch 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. Received on Fri Aug 18 00:54:23 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Fri 18 Aug 2006 - 02:22:14 EST.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.