[R] fitting truncated normal distribution

Sundar Dorai-Raj sundar.dorai-raj at pdf.com
Fri Aug 18 14:13:11 CEST 2006


Hi, Markus,

One other suggestion is to add the "lower" argument to fitdistr:

fitdistr(x, dtnorm0, start = list(mean = 0, sd = 1), lower = 0)

where dtnorm0 is defined as before. This indicates to fitdistr that the 
optimization should be constrained. See ?optim for details.

--sundar

Schweitzer, Markus wrote:
>  Thank you Sundar,
> 
> Yes, always integers. By demand data I meant the amount of ordered
> products in a certain period. Therefore, x is a vector of periods (i.e.
> Weeks in a year)
> 
> In my example we could see an article, that has only been ordered in two
> weeks within one year.
> All the zeros show, that nobody has ordered the items in these periods.
> (First half of the year/first 24 weeks)
> 
> Since the orders cannot be negative, some literature recommended to use
> a truncated normal distribution (Poisson and negative binomial are also
> recommended).
> 
> My "x" is just a sample out of the dataset. There might be other time
> series with better attributes for a truncated normal distribution.
> 
> My problem is simply, that I only get an error message when I use
> fitdistr.
> 
> 
>>>>"Error in "[<-"(`*tmp*`, x >= lower & x <= upper, value = numeric(0))
> 
> nothing to replace.
> 
> I hope, there is a way fitdistr can also compute "difficult" data.
> 
> Best regards, markus
> 
> 
> 
> -----Original Message-----
> From: Sundar Dorai-Raj [mailto:sundar.dorai-raj at pdf.com] 
> Sent: Donnerstag, 17. August 2006 16:47
> To: Schweitzer, Markus
> Cc: r-help at stat.math.ethz.ch
> Subject: Re: [R] fitting truncated normal distribution
> 
> 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 at 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 at 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 at 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 at 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 at 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 at 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.



More information about the R-help mailing list