R help - Oct 2012 - Help with R Fitting an inverse Gamma

Dear all,

I am new in R and would like to ask for someone's help in understanding
where I go wrong with the following code:

rm(list=ls())
# Required packages
library(MCMCpack)

# Simulated data
set.seed(1)
data = rinvgamma(n=250, shape = 5, scale = 2) + 2

hist(data)

# log-likelihood
ll = function(par){
if(par[1]>0 & par[2]>0 & par[3]<min(data)) return(
-sum(log(dinvgamma(data-
par[3],par[1],par[2]))) )
else return(Inf)
}

# MLE
mle = optim(c(5,2,2),ll)
params = mle$par

# Fit
hist(data,probability=T,ylim=c(0,2.5))
points(seq(2,4.5,0.001),dinvgamma(seq(2,4.5,0.001)-params[3],params[1],params[2]),type=&quot;l&quot;,col=&quot;red&quot;)

This code fits an Inverse Gamma distribution to the randomly generated data
and plots the associated histogram. My problem is that when I try to import
a dataset saved locally it gives me an error. Here is the modification of
the code: 



# Required packages
library(MCMCpack)

data=my.csv.data$V1
hist(data)

# log-likelihood
ll = function(par){
if(par[1]>0 & par[2]>0 & par[3]<min(data)) return(
-sum(log(dinvgamma(data-
par[3],par[1],par[2]))) )
else return(Inf)
}

# MLE
mle = optim(c(5,2,2),ll)
params = mle$par

# Fit
hist(data,probability=T,ylim=c(0,2.5))
points(seq(2,4.5,0.001),dinvgamma(seq(2,4.5,0.001)-params[3],params[1],params[2]),type=&quot;l&quot;,col=&quot;red&quot;)



I have attached the file which I import directly from the menu, named:
excel250.csv. When I run the code it gives me the following error:

Error in optim(c(5, 2, 2), ll) : 
  function cannot be evaluated at initial parameters

I would appreciate your help.

Thank you.

Konstantinos

&lt;nabble_a href=&quot;excel250.csv&quot;>excel250.csv 



--
View this message in context:
http://r.789695.n4.nabble.com/Help-with-R-Fitting-an-inverse-Gamma-tp4644984.html
Sent from the R help mailing list archive at Nabble.com.

On Thu, Oct 4, 2012 at 11:42 AM, kmammasis <mammasis82 at hotmail.com>
wrote:
> Dear all,
>
> I am new in R and would like to ask for someone's help in understanding
> where I go wrong with the following code:
>
# > rm(list=ls())
#

Please don't include this line into your postings: it's quite rude to
clobber the data of folks trying to help you.

> # Required packages
> library(MCMCpack)
>
> # Simulated data
> set.seed(1)
> data = rinvgamma(n=250, shape = 5, scale = 2) + 2
>
# Bad name for a data variable: see
library(fortunes)
fortune("dog")
> hist(data)
>
> # log-likelihood
> ll = function(par){
> if(par[1]>0 & par[2]>0 & par[3]<min(data)) return(
-sum(log(dinvgamma(data-
> par[3],par[1],par[2]))) )
> else return(Inf)
> }
>
> # MLE
> mle = optim(c(5,2,2),ll)
> params = mle$par
>
> # Fit
> hist(data,probability=T,ylim=c(0,2.5))
>
points(seq(2,4.5,0.001),dinvgamma(seq(2,4.5,0.001)-params[3],params[1],params[2]),type="l",col="red")
>
> This code fits an Inverse Gamma distribution to the randomly generated data
> and plots the associated histogram. My problem is that when I try to import
> a dataset saved locally it gives me an error. Here is the modification of
> the code:
>
>
>
> # Required packages
> library(MCMCpack)
>
> data=my.csv.data$V1
> hist(data)
>
> # log-likelihood
> ll = function(par){
> if(par[1]>0 & par[2]>0 & par[3]<min(data)) return(
-sum(log(dinvgamma(data-
> par[3],par[1],par[2]))) )
> else return(Inf)
> }
>
> # MLE
> mle = optim(c(5,2,2),ll)
> params = mle$par
>
> # Fit
> hist(data,probability=T,ylim=c(0,2.5))
>
points(seq(2,4.5,0.001),dinvgamma(seq(2,4.5,0.001)-params[3],params[1],params[2]),type="l",col="red")
>
>
>
> I have attached the file which I import directly from the menu, named:
> excel250.csv. When I run the code it gives me the following error:
No you haven't. The R-help servers scrub most attachments. Use dput()
to make a nice reproducible example:

http://stackoverflow.com/questions/5963269/how-to-make-a-great-r-reproducible-example

Cheers,
M
>
> Error in optim(c(5, 2, 2), ll) :
>   function cannot be evaluated at initial parameters
Totally random guess: your function isn't being evaluated properly
either because i) you used the name data() which is also a function
and R is having trouble finding your data; ii) your data is such that
dinvgamma() gives bad value at the initial point c(5,2,2)

Cheers,
Michael
>
> I would appreciate your help.
>
> Thank you.
>
> Konstantinos
>
> <nabble_a href="excel250.csv">excel250.csv
>
>
>
> --
> View this message in context:
http://r.789695.n4.nabble.com/Help-with-R-Fitting-an-inverse-Gamma-tp4644984.html
> Sent from the R help mailing list archive at Nabble.com.
>
> ______________________________________________
> R-help at r-project.org 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 and thanks for your replies. I have changed the name of the input data to
something other than "data" to avoid confusion with the function
data(),
however the problem still persists. I get the same error. I have found it
easier to upload the file in the following location:

http://www.sendspace.com/file/euv6qy
<http://www.sendspace.com/file/euv6qy>

After importing the csv file in R I use the following code:

# Required packages
library(MCMCpack)

V=my.csv.data$V1
hist(V)

# log-likelihood
ll = function(par){
if(par[1]>0 & par[2]>0 & par[3]<min(V)) return(
-sum(log(dinvgamma(V-
par[3],par[1],par[2]))) )
else return(Inf)
}

# MLE
mle = optim(c(0.05,.022,.0sessionInfo()2),ll)
params = mle$par

# Fit
hist(Vdput()
,probability=T,ylim=c(0,2.5))
points(seq(2,4.5,0.001),dinvgamma(seq(2,4.5,0.001)-params[3],params[1],params[2]),type="l",col="red")

Thank you once more.

BR,

Konstantinos





--
View this message in context:
http://r.789695.n4.nabble.com/Help-with-R-Fitting-an-inverse-Gamma-tp4644984p4645128.html
Sent from the R help mailing list archive at Nabble.com.