R help - Aug 2006 - fitting truncated normal distribution

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)

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

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]]

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.