Why to optimize max log probability instead of probability In most machine learning tasks where you can formulate some probability $p$ which should be maximised, we would actually optimize the log probability $\log p$ instead of the probability for some parameters $\theta$. E.g. in maximum likelihood training, it's usually the log-likelihood. When doing this with some gradient method, this involves a factor:
$$ \frac{\partial \log p}{\partial \theta} = \frac{1}{p} \cdot \frac{\partial p}{\partial \theta} $$
See here or here for some examples.
Of course, the optimization is equivalent, but the gradient will be different, so any gradient-based method will behave different (esp. stochastic gradient methods).
Is there any justification that the $\log p$ gradient works better than the $p$ gradient?
 A: Underflow
The computer uses a limited digit floating point representation of fractions, multiplying so many probabilities is guaranteed to be very very close to zero.
With $log$, we don't have this issue.
A: *

*The logarithm of the probability of multiple joint probabilities simplifies to the sum of the logarithms of the individual probabilities (and the sum rule is easier than the product rule for differentiation)
$\log \left(\prod_i P(x_i)\right) = \sum_i \log \left( P(x_i)\right)$

*The logarithm of a member of the family of exponential probability distributions (which includes the ubiquitous normal) is polynomial in the parameters (i.e. max-likelihood reduces to least-squares for normal distributions)
$\log\left(\exp\left(-\frac{1}{2}x^2\right)\right) = -\frac{1}{2}x^2$

*The latter form is both more numerically stable and symbolically easier to differentiate than the former.

*Last but not least, the logarithm is a monotonic transformation that preserves the locations of the extrema (in particular, the estimated parameters in max-likelihood are identical for the original and the log-transformed formulation)
A: By using the $\ln p$ we increase the dynamic range of the optimization algorithm. The $p$ in applications is usually a product of functions. For instance, in maximum likelihood estimation it's the product of the form $L(x|\theta)=\Pi_{i=1}^n f(x_i|\theta)$, where $f(.)$ is the density function, which can be greater or less than 1, btw.
So, when $n$ is very large, i.e. large sample, your likelihood function $L(.)$ is usually far away from 1: it's either very small or very large, because it's a power function $L\sim f(.)^n$. 
By taking a log we simply improve the dynamic range of any optimization algorithm, allowing it to work with extremely large or small values in the same way.
A: It is much easier to take a derivative of sum of logarithms than to take a derivative of product, that contains, say, 100 multipliers.
A: As a general rule, the most basic and easy optimization problem is to optimize a quadratic function. You can easily find the optimum of such a function no matter where you start. How this manifests depends on the specific method but the closer your function to a quadratic, the better.
As noted by TemplateRex, in a wide variety of problems, the probabilities that go into calculating the likelihood function come from the normal distribution, or are approximated by it. So if you work on the log, you get a nice quadratic function. Whereas if you work on the probabilities, you have a function that


*

*Is not convex (the bane of optimization algorithms everywhere)

*Crosses multiple scales rapidly, and therefore has a very narrow range where the function values are indicative of where to direct your search.


Which function would you rather optimize, this, or this?
(This was actually an easy one; in practical applications your search can start so far off the optimum that the function values and gradients, even if you were able to compute them numerically, will be indistinguishable from 0 and useless for the purposes of the optimization algorithm. But transforming to a quadratic function makes this a piece of cake.)
Note that this is completely consistent with the numerical stability issues already mentioned. The reason log scale is required to work with this function, is exactly the same reason that the log probability is much better behaved (for optimization and other purposes) than the original.
You could also approach this another way. Even if there was no advantage to the log (which there is) - we're gonna use the log scale anyway for derivations and calculation, so what reason is there to apply the exp transformation just for computing the gradient? We may as well remain consistent with the log.
A: Gradient methods generally work better optimizing $\log p(x)$ than $p(x)$ because the gradient of $\log p(x)$ is generally more well-scaled. That is, it has a size that consistently and helpfully reflects the objective function's geometry, making it easier to select an appropriate step size and get to the optimum in fewer steps.
To see what I mean, compare the gradient optimization process for $p(x) = \exp(-x^2)$ and $f(x) = \log p(x) = -x^2$. At any point $x$, the gradient of $f(x)$ is $$f'(x) = -2x.$$ If we multiply that by $1/2$, we get the exact step size needed to get to the global optimum at the origin, no matter what $x$ is. This means that we don't have to work too hard to get a good step size (or "learning rate" in ML jargon). No matter where our initial point is, we just set our step to half the gradient and we'll be at the origin in one step. And if we don't know the exact factor that is needed, we can just pick a step size around 1, do a bit of line search, and we'll find a great step size very quickly, one that works well no matter where $x$ is. This property is robust to translation and scaling of $f(x)$. While scaling $f(x)$ will cause the optimal step scaling to differ from 1/2, at least the step scaling will be the same no matter what $x$ is, so we only have to find one parameter to get an efficient gradient-based optimization scheme.
In contrast, the gradient of $p(x)$ has very poor global properties for optimization. We have $$p'(x) = f'(x) p(x)= -2x \exp(-x^2).$$ This multiplies the perfectly nice, well-behaved gradient $-2x$ with a factor $\exp(-x^2)$ which decays (faster than) exponentially as $x$ increases. At $x = 5$, we already have $\exp(-x^2) = 1.4 \cdot 10^{-11}$, so a step along the gradient vector is about $10^{-11}$ times too small. To get a reasonable step size toward the optimum, we'd have to scale the gradient by the reciprocal of that, an enormous constant $\sim 10^{11}$. Such a badly-scaled gradient is worse than useless for optimization purposes - we'd be better off just attempting a unit step in the uphill direction than setting our step by scaling against $p'(x)$! (In many variables $p'(x)$ becomes a bit more useful since we at least get directional information from the gradient, but the scaling issue remains.)
In general there is no guarantee that $\log p(x)$ will have such great gradient scaling properties as this toy example, especially when we have more than one variable. However, for pretty much any nontrivial problem, $\log p(x)$ is going to be way, way better than $p(x)$. This is because the likelihood is a big product with a bunch of terms, and the log turns that product into a sum, as noted in several other answers. Provided the terms in the likelihood are well-behaved from an optimization standpoint, their log is generally well-behaved, and the sum of well-behaved functions is well-behaved. By well-behaved I mean $f''(x)$ doesn't change too much or too rapidly, leading to a nearly quadratic function that is easy to optimize by gradient methods. The sum of a derivative is the derivative of the sum, no matter what the derivative's order, which helps to ensure that that big pile of sum terms has a very reasonable second derivative!
A: Some nice answers have been given already. But I encountered recently a new one:
Often, you are given a huge training data set $\mathcal{X}$, and you define some probabilistic model $p(x|\theta)$, and you want to maximize the likelihood for $x \in \mathcal{X}$. They are assumed to be independent, i.e. you have
$$
p(\mathcal{X}|\theta) = \prod_{x\in\mathcal{X}} p(x|\theta) .
$$
Now, you often do some sort of stochastic (mini-batch) gradient-based training, i.e. in each step, for your loss $L$, you optimize $L(\mathcal{X'}|\theta)$ for $\mathcal{X'} \subset \mathcal{X}$, i.e.
$$
\theta' := \theta - \frac{\partial \sum_{x\in\mathcal{X'}} L(x|\theta)}{\partial \theta} .
$$
Now, these stochastic steps are accumulated additively. Because of that, you want the property that in general
$$
L(\mathcal{X}|\theta) = \sum_{x\in\mathcal{X}} L(x|\theta) .
$$
This is the case for
$$
L(x|\theta) = -\log p(x|\theta) .
$$
