[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