[R] Help on R commands
Ben Bolker
ben at zoo.ufl.edu
Fri Dec 20 13:13:02 CET 2002
You're not too far off.
You need to do something like
objfun <- function(p) {
mu <- p[1]
sigma <- p[2]
-sum(dnorm(x,mean=mu,sd=sigma,log=TRUE))
}
optim(par=c(0,1),fn=objfun) # check ?optim for correct argument names
Some points:
- the main thing is that you have to define mu and sigma within your
objective function
- you're better off with the log-likelihood than the likelihood; by
convention people usually try to minimize the negative log-likelihood
rather than maximizing the likelihood
- "nlmin" is an S-PLUS function. The analogues in R are nlm and optim.
- your setting negative values to zero complicates your problem
considerably. Can you use a gamma model instead?
- the ML estimate of normal data (given that you haven't messed with
negative values) is mu=mean(x), sd=sd(x)*sqrt(length(x)/(length(x)-1)) =
sqrt(mean((x-mean(x))^2))
A lot of these questions are pretty basic; you really need to try to get
some help off-list. You may not get too much more help here unless you
explain what your situation is (are you doing this for a class,
independent project, etc.) and give some indication that you've exhausted
local resources ...
cheers,
Ben Bolker
On Wed, 18 Dec 2002 Jameshutton25 at aol.com wrote:
> Dear mailing list,
>
> I desperately need help on making a small program that is trying to find the
> likelihood of a distribution. Anyone that has any ideas please feel free to
> suggest them.
>
> Ok this is what I have done so far:
>
> I wanted 20 random numbers that were normally distributed, and I did this by
> typing x<-rnorm(20).
> I then wanted to change any negative values in the data set to zero and I
> did this by x[x<0]<-0.
> These numbers have come from the normal density
> 1/sigma(2pi)1/2exp-{(x-mew)2/2sigma2}, what I want to do now and having
> trouble with is that for each of these results (which can be substituted back
> into the x in the eqation above) is to multiply all 20 results (in the above
> equation form) together to form the likelihood. However the critical problem
> that I am experiencing is that for each case I do not know the mew and sigma
> eventhough that i know tthe x value. The reason why i dont know the vales of
> sigma and mew is because after I have formed the program i want to use nlmin(
> , ) to basically maximse the the likelihood so i find the values for mew
> and sigma, this is what i am aiming finally to do.
> My pathetic effort so far is:
>
> prod(dnorm(x,mean=mu,sd=sigma)
>
> However I know that this doesn't incorporate the fact that mew and sigma are
> not known as when i input this it says that mu is not recognised but I dont
> know how to make mu and sigma different to each other.
>
> If anyone has any ideas please feel free to suggest them as I will basically
> try anything.
>
> Yours Faithfully,
>
> James Hutton
>
> [[alternate HTML version deleted]]
>
> ______________________________________________
> R-help at stat.math.ethz.ch mailing list
> http://www.stat.math.ethz.ch/mailman/listinfo/r-help
>
--
318 Carr Hall bolker at zoo.ufl.edu
Zoology Department, University of Florida http://www.zoo.ufl.edu/bolker
Box 118525 (ph) 352-392-5697
Gainesville, FL 32611-8525 (fax) 352-392-3704
More information about the R-help
mailing list