I need to optimize parameters of a Gaussian shaped function to best fit my data points using maximum likelihood.  I first make an initial estimate for the parameters (mu, sigma and peak) and apply a gradient descent (fmincon in Matlab) algorithm to optimize these parameters.   The cost function is the sum of the negative log values of probabilities of observing each data point given the model parameters and the model noise parameter.   sum(-log( normpdf( Y - model( X, params) ,0,error_sigma) ));   where, Y: observations  normpdf : normal probability densitiy function.  model:  a function that takes the x values and parameters to be optimized  params : parameters of the model  error_sigma : noise of the data  I optimize the model parameters and the noise parameter simultaneously to obtain highest -loglikelihood values.  I am puzzled by obtaining negative likelihood values. I would have estimated that -log of probabilities should always evaluate in positive values. However, the matlab function normpdf, especially for small spread values, returns probability values which are bigger than 1.  This is because the integral of the normpdf function does nearly never sums up to 1. Therefore I am thinking that I should be using the normal cumulative function for optimization, instead of normpdf in order to get valid loglikelihood values. However, the parameters never converges as good as the normpdf case, when I use normcdf.  I hope I am clear.
470
On the , there is this line:      $p(\mathrm{height}|\mathrm{male}) = 1.5789$ (A probability distribution over 1 is OK. It is the area under the bell curve that is equal to 1.)    How can a value $&gt;1$ be OK? I thought all probability values were expressed in the range $0 \leq p \leq 1$. Furthermore, given that it is possible to have such a value, how is that value obtained in the example shown on the page?
1419
There is a classic problem:  Suppose that $X_1,\ldots,X_n$ form an i.i.d. sample from a uniform distribution on the interval $(0,\theta)$, where $\theta&gt;0$ is unknown. I would like to find the MLE of $\theta$.  The pdf of each observation will have the form: $$ f(x\mid\theta) =  \begin{cases} 1/\theta\quad&amp;\text{for }\, 0\leq x\leq \theta\\ 0 &amp;\text{otherwise}. \end{cases} $$ The likelihood function therefore has the form: $$ L(\theta) = \begin{cases} 1/\theta^n \quad&amp;\text{for }\; 0\leq x_i \leq \theta\;\; \text{for all }i,\\ 0 &amp;\text{otherwise}. \end{cases} $$ The general solution is usually that the MLE of theta must be a value of $\theta$ for which $\theta \geq x_i$ and which maximizes $1/\theta^n$ among all such values.   The reasoning is that since $1/\theta^n$ is a decreasing function of $\theta$, the estimate will be the smallest possible value of $\theta$ such that $\theta\geq x_i$.  Therefore, the mle of $\theta$, $\hat{\theta}$, is $\max(X_1,\ldots,X_n)$.  Here, I do not understand why we cannot just differentiate the likelihood function with respect to theta and then set it equal to $0$?  Thanks!
102503
A participant in one experiment needs to decide whether a flash and a sound are simultaneous or not for many possible asynchronies between the flash and the sound (x in seconds). For each asynchrony, the flash and sound are presented 100 times. In the below graph the proportion of 'simultaneous' responses (y) is plotted as a function of the asynchrony    I want to fit a 3-parameter Gaussian distribution to these data. I used least squares (below). My question is whether least squares is the most standard/orthodox form to proceed for this kind of fitting. It makes sense to use maximum likelihood estimation here? How can be done?     x&lt;-seq(-.3,.3,.1)   yes&lt;-c(10,45,90,100,60,10,5) #simultaneous responses   y&lt;-yes/100    leastSquares&lt;-function(p) sum((y-p[1]*dnorm(x,p[2],p[3]))^2)   p&lt;-optim(c(.001,0.1,.210),leastSquares)$par    xseq&lt;-seq(-.3,.3,.01)   yseq&lt;-p[1]*dnorm(xseq,p[2],p[3])    plot(x,y)   lines(xseq,yseq)
74277
Hi can the same be shown to obtain shape and scale parameter for modified maximum likelihood method
58721
I have a set of observation, let's call it $X$ and would like to fit a cdf to it. $X$ has a distribution which is roughly approximable with the normal distribution. This CDF should correspond to a continuous distribution function.  So far I've used a parametric approach by estimating mean and standard deviation and using a normal cdf but I would like to know what other options are available and how to use them.  How does the set of available option change if I require the cdf to be a smooth curve?
95513
Let $X_1, X_2, ..X_n$ represent a random sample from this pdf:  $$f(x|\theta)= \frac{3x^2}{\theta^3}, 0\leq x\leq \theta$$ (with 0 elsewhere)  Could someone explain how I would find the maximum likelihood estimator (MLE) of this pdf? I know I have to get the likelihood function L, then the log likelihood function, then differentiate, but I am getting lost in the calculations. Thanks!
118652
Let $X_1,\ldots,X_n$ be a sample from the distribution whose probability density   function is:       $$f(x)=\begin{cases}\frac13e^{-\frac13(x-\theta)}&amp;, x\ge\theta \\ 0&amp;,x&lt;\theta\end{cases}$$      Find the maximum likelihood estimator of $\theta$.   So the PDF is equal to:  $$f(x\mid \theta) = \frac{1}{3}e^{-\frac{1}{3}(x-\theta)}1_{\{x\geq \theta\}}$$  and thereby the likelihood will be:   $$\mathcal{L} = \prod_{i=1}^n f(x\mid \theta)= \prod_{i=1}^n  \frac{1}{3}e^{-\frac{1}{3}(x-\theta)}1_{\{x\geq \theta\}}$$
34102
I want to define a Gaussian distribution function and plot it in python using the mode and inflection points parameter values instead of using the mean and standard deviation. For example, I have mode=110 and two points : (40, 160) for asymmetrical points. Or mode=100 and two points = (50,150) for symmetrical points. I tried with this code below, but I don't how I can add the second point to the formula which is returned by the function? Since I have two points, and I don't know how to define a function that depends on two points and mode. Thanks in advance for your help! import numpy as np matplotlib.pyplot as plt    def gaussian(x, mode, inf_point):     return 1/(np.sqrt(2*np.pi)*inf_point)*np.exp(-np.power((x - mode)/inf_point, 2)/2) x = np.linspace(0,256) plt.plot(x, gaussian(x, mode, inf_point))
123839
I am trying to get the maximum likelihood estimate for the parameter $p$. The distribution is the following:   $$ f(x\mid p) = \begin{cases} \frac{p}{x^2} &amp;\text{for} \ p\leq x &lt; \infty  \ \\ 0 &amp;\text{if not} \end{cases} $$  The sample has size $n$.  The problem is, when I try to estimate it by the procedure I know, I would have to estimate the likelihood function and obtain the derivative of the log-likelihood. We'd have:   $$ L(p; x) = \frac{p^n}{\prod_{i=1}^{n} x^2_i}$$  $$ \ln L(p;x) = n \ln(p) - \sum_{i=1}^{n} \ln(x_i^2) $$   For the derivative:  $$ l'(p;x) = \frac{n}{p} = 0$$   And I am stuck because it has no solution for $p$. How do I evaluate this?  Thanks!   EDIT:  So in this case I can use the indicator variable to write:   $$ L(p;x) = \frac{p^n}{\prod_{i=1}^{n}x_i^2} I_{(x_i \geq p)}$$   for $i = 1,2, \ldots n$ in the indicator variable. So the "closest" non-null value of $x \in X$ to $p$ is $min(X_1, \ldots X_n)$. Is that the point?
47871
Suppose we have a random variable $X \sim f(x|\theta)$. If $\theta_0$ were the true parameter, the the likelihood function should be maximized and the derivative equal to zero. This is the basic principle behind the maximum likelihood estimator.   As I understand it, Fisher information is defined as   $$I(\theta) = \Bbb E \Bigg[\left(\frac{\partial}{\partial \theta}f(X|\theta)\right)^2\Bigg ]$$  Thus, if $\theta_0$ is the true parameter, $I(\theta) = 0$. But if it $\theta_0$ is not the true parameter, then we will have a larger amount of Fisher information.   my questions    Does Fisher information measure the "error" of a given MLE? In other words, doesn't the existence of positive Fisher information imply my MLE can't be ideal?  How does this definition of "information" differ from that used by Shannon? Why do we call it information?
105658
I'm reading about maximum likelihood .  In the last paragraph of the first page, it says:    Why does the value of $p$ that maximizes $\log L(p;3)$ is the same $p$ that maximizes $L(p;3)$. The fact that it mentions that log is an increasing function does not help me understand it at all. Does this method work for any increasing function?  I've covered single and two variable optimization but have not come across this fact until now. I'd like to understand why it must be true and read more about its explanation. So, I would also appreciate it if someone can provide some links to it (I don't know what term to search).
123437
Suppose $X_1, X_2, \dots, X_n$ are iid Bernoulli(p) random variables. How do you find the restricted maximum likelihood for p where $0&lt;p&lt;0.5$?  My work so far:  Write out the likelihood: $L=p^{\sum x} (1-p)^{n-\sum x}$  Take the log: $\ln(L) = \sum x\ln(p)+(n-\sum x)\ln(1-p)$  Take the derivative and set equal to zero: $\frac{\sum x-np}{p(1-p)}=0$.  Now I'm stuck. I don't know how to incorporate the restrictions for $p$.
23509
Suppose a differentiable, convex function $F(x)$ exists. Then $b = a - \gamma\nabla F(a)$ implies that $F(b) \leq F(a)$ given $\gamma$ is chosen properly. The goal is to find the optimal $\gamma$ at each step. In my book, in order to do this, one should minimize  $G(\gamma)=$ $F(x-\gamma\nabla F(x))$ for $\gamma$. It also says that it should be minimized via a line search. Why can't this function be minimized by simple calculus?  I am unsure what it means to perform a line search on this function. I am also unsure why we want to minimize this function to find the optimal step size.
4299
So in a normal distribution, we have two parameters: mean $\mu$ and variance $\sigma^2$. In the book Pattern Recognition and Machine Learning, there suddenly appears a hyperparameter $\lambda$ in the  regularization terms of the error function.  What are hyperparameters? Why are they named as such? And how are they intuitively different from parameters in general?
121545
In MLE for continuous case, my course notes define the likelihood function to be:  $$ L(\theta) = L(\theta;y) = \prod_{i=1}^n f(y_i;\theta) $$  Where $f$ is the joint pdf of $y_i$ given $\theta$.  I thought normally it doesn't make sense to evaluate a pdf at a point (eg. $f(X=2) = 0$ for some continuous rv X), and that it only make sense to evaluate for the random variable to be $&gt;$ or $&lt;$ to a value.  But in the case of MLE, why is it OK to evaluate a joint pdf at a point?
41895
I try to reproduce with optim the results from a simple linear regression fitted with glm or even nls R functions. The parameters estimates are the same but the residual variance estimate and the standard errors of the other parameters are not the same particularly when the sample size is low. I suppose that this is due differences in the way the residual standard error is calculated between Maximum Likelihood and Least square approaches (dividing by n or by n-k+1 see bellow in the example). I understand from my readings on the web that optimization is not a simple task but I was wondering if it would be possible to reproduce in a simple way the standard error estimates from glm while using optim.   Simulate a small dataset  set.seed(1) n = 4 # very small sample size ! b0 &lt;- 5 b1 &lt;- 2 sigma &lt;- 5 x &lt;- runif(n, 1, 100) y =  b0 + b1*x + rnorm(n, 0, sigma)    Estimate with optim  negLL &lt;- function(beta, y, x) {     b0 &lt;- beta[1]     b1 &lt;- beta[2]     sigma &lt;- beta[3]     yhat &lt;- b0 + b1*x     likelihood &lt;- dnorm(y, yhat, sigma)     return(-sum(log(likelihood))) }  res &lt;- optim(starting.values, negLL, y = y, x = x, hessian=TRUE) estimates &lt;- res$par     # Parameters estimates se &lt;- sqrt(diag(solve(res$hessian))) # Standard errors of the estimates cbind(estimates,se)       &gt; cbind(estimates,se)       estimates         se b0     9.016513 5.70999880 b1     1.931119 0.09731153 sigma  4.717216 1.66753138   Comparison with glm and nls  &gt; m &lt;- glm(y ~ x) &gt; summary(m)$coefficients             Estimate Std. Error   t value    Pr(&gt;|t|) (Intercept) 9.016113  8.0759837  1.116411 0.380380963 x           1.931130  0.1376334 14.030973 0.005041162 &gt; sqrt(summary(m)$dispersion) # residuals standard error [1] 6.671833 &gt;  &gt; summary(nls( y ~ b0 + b1*x, start=list(b0 = 5, b1= 2)))  Formula: y ~ b0 + b1 * x  Parameters:    Estimate Std. Error t value Pr(&gt;|t|)    b0   9.0161     8.0760   1.116  0.38038    b1   1.9311     0.1376  14.031  0.00504 ** --- Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1  Residual standard error: 6.672 on 2 degrees of freedom   I can reproduce the different residual standard error estimates like this :   &gt; # optim / Maximum Likelihood estimate &gt; sqrt(sum(resid(m)^2)/n) [1] 4.717698 &gt;  &gt; # Least squares estimate (glm and nls estimates) &gt; k &lt;- 3 # number of parameters &gt; sqrt(sum(resid(m)^2)/(n-k+1)) [1] 6.671833
77818
Let $A,B,C,D$ be independent normal random variables distributed according to $\mathcal{N}(0,\sigma^2)$ and $E = A+B+C+D$.  The goal is to provide an optimal estimate of $\sigma^2$ given samples of $A,B,E$ but not $C,D$.  In this case optimal refers to minimising the root-mean square deviation of the estimate $\hat{\sigma}^2$:  $$RMSD(\hat{\sigma}^2) = \sqrt{\mathbb{E}\left[\left(\hat{\sigma}^2-\sigma^2\right)^2\right]}$$  I can think of a simple approach to estimate $\sigma^2$ given samples of $A,B$ which is to simply compute the sample variance. Also I know that the variance of $E$ is $4\sigma^2$ so I could compute the sample variance of $E$ and divide it by 4 to get an estimate of $\sigma^2$. However it should be possible to combine the information from both $E$ and $A,B$ to improve estimation $\sigma^2$. Could you outline the approach to use both sets of information?